QGIS API Documentation  3.6.0-Noosa (5873452)
qgsgeometrysnapper.cpp
Go to the documentation of this file.
1 /***************************************************************************
2  * qgsgeometrysnapper.cpp *
3  * ------------------- *
4  * copyright : (C) 2014 by Sandro Mani / Sourcepole AG *
5  * email : [email protected] *
6  ***************************************************************************/
7 
8 /***************************************************************************
9  * *
10  * This program is free software; you can redistribute it and/or modify *
11  * it under the terms of the GNU General Public License as published by *
12  * the Free Software Foundation; either version 2 of the License, or *
13  * (at your option) any later version. *
14  * *
15  ***************************************************************************/
16 
17 #include <QtConcurrentMap>
18 #include "qgsfeatureiterator.h"
19 #include "qgsgeometry.h"
20 #include "qgsvectorlayer.h"
21 #include "qgsgeometrysnapper.h"
22 #include "qgsvectordataprovider.h"
23 #include "qgsgeometryutils.h"
24 #include "qgsmapsettings.h"
25 #include "qgssurface.h"
26 #include "qgsmultisurface.h"
27 #include "qgscurve.h"
28 
30 
31 QgsSnapIndex::PointSnapItem::PointSnapItem( const QgsSnapIndex::CoordIdx *_idx, bool isEndPoint )
32  : SnapItem( isEndPoint ? QgsSnapIndex::SnapEndPoint : QgsSnapIndex::SnapPoint )
33  , idx( _idx )
34 {}
35 
36 QgsPoint QgsSnapIndex::PointSnapItem::getSnapPoint( const QgsPoint &/*p*/ ) const
37 {
38  return idx->point();
39 }
40 
41 QgsSnapIndex::SegmentSnapItem::SegmentSnapItem( const QgsSnapIndex::CoordIdx *_idxFrom, const QgsSnapIndex::CoordIdx *_idxTo )
42  : SnapItem( QgsSnapIndex::SnapSegment )
43  , idxFrom( _idxFrom )
44  , idxTo( _idxTo )
45 {}
46 
47 QgsPoint QgsSnapIndex::SegmentSnapItem::getSnapPoint( const QgsPoint &p ) const
48 {
49  return QgsGeometryUtils::projectPointOnSegment( p, idxFrom->point(), idxTo->point() );
50 }
51 
52 bool QgsSnapIndex::SegmentSnapItem::getIntersection( const QgsPoint &p1, const QgsPoint &p2, QgsPoint &inter ) const
53 {
54  const QgsPoint &q1 = idxFrom->point(), & q2 = idxTo->point();
55  QgsVector v( p2.x() - p1.x(), p2.y() - p1.y() );
56  QgsVector w( q2.x() - q1.x(), q2.y() - q1.y() );
57  double vl = v.length();
58  double wl = w.length();
59 
60  if ( qgsDoubleNear( vl, 0, 0.000000000001 ) || qgsDoubleNear( wl, 0, 0.000000000001 ) )
61  {
62  return false;
63  }
64  v = v / vl;
65  w = w / wl;
66 
67  double d = v.y() * w.x() - v.x() * w.y();
68 
69  if ( d == 0 )
70  return false;
71 
72  double dx = q1.x() - p1.x();
73  double dy = q1.y() - p1.y();
74  double k = ( dy * w.x() - dx * w.y() ) / d;
75 
76  inter = QgsPoint( p1.x() + v.x() * k, p1.y() + v.y() * k );
77 
78  double lambdav = QgsVector( inter.x() - p1.x(), inter.y() - p1.y() ) * v;
79  if ( lambdav < 0. + 1E-8 || lambdav > vl - 1E-8 )
80  return false;
81 
82  double lambdaw = QgsVector( inter.x() - q1.x(), inter.y() - q1.y() ) * w;
83  return !( lambdaw < 0. + 1E-8 || lambdaw >= wl - 1E-8 );
84 }
85 
86 bool QgsSnapIndex::SegmentSnapItem::getProjection( const QgsPoint &p, QgsPoint &pProj )
87 {
88  const QgsPoint &s1 = idxFrom->point();
89  const QgsPoint &s2 = idxTo->point();
90  double nx = s2.y() - s1.y();
91  double ny = -( s2.x() - s1.x() );
92  double t = ( p.x() * ny - p.y() * nx - s1.x() * ny + s1.y() * nx ) / ( ( s2.x() - s1.x() ) * ny - ( s2.y() - s1.y() ) * nx );
93  if ( t < 0. || t > 1. )
94  {
95  return false;
96  }
97  pProj = QgsPoint( s1.x() + ( s2.x() - s1.x() ) * t, s1.y() + ( s2.y() - s1.y() ) * t );
98  return true;
99 }
100 
102 
103 class Raytracer
104 {
105  // Raytrace on an integer, unit-width 2D grid with floating point coordinates
106  // See http://playtechs.blogspot.ch/2007/03/raytracing-on-grid.html
107  public:
108  Raytracer( float x0, float y0, float x1, float y1 )
109  : m_dx( std::fabs( x1 - x0 ) )
110  , m_dy( std::fabs( y1 - y0 ) )
111  , m_x( std::floor( x0 ) )
112  , m_y( std::floor( y0 ) )
113  , m_n( 1 )
114  {
115  if ( m_dx == 0. )
116  {
117  m_xInc = 0.;
118  m_error = std::numeric_limits<float>::infinity();
119  }
120  else if ( x1 > x0 )
121  {
122  m_xInc = 1;
123  m_n += int( std::floor( x1 ) ) - m_x;
124  m_error = ( std::floor( x0 ) + 1 - x0 ) * m_dy;
125  }
126  else
127  {
128  m_xInc = -1;
129  m_n += m_x - int( std::floor( x1 ) );
130  m_error = ( x0 - std::floor( x0 ) ) * m_dy;
131  }
132  if ( m_dy == 0. )
133  {
134  m_yInc = 0.;
135  m_error = -std::numeric_limits<float>::infinity();
136  }
137  else if ( y1 > y0 )
138  {
139  m_yInc = 1;
140  m_n += int( std::floor( y1 ) ) - m_y;
141  m_error -= ( std::floor( y0 ) + 1 - y0 ) * m_dx;
142  }
143  else
144  {
145  m_yInc = -1;
146  m_n += m_y - int( std::floor( y1 ) );
147  m_error -= ( y0 - std::floor( y0 ) ) * m_dx;
148  }
149  }
150  int curCol() const { return m_x; }
151  int curRow() const { return m_y; }
152  void next()
153  {
154  if ( m_error > 0 )
155  {
156  m_y += m_yInc;
157  m_error -= m_dx;
158  }
159  else if ( m_error < 0 )
160  {
161  m_x += m_xInc;
162  m_error += m_dy;
163  }
164  else
165  {
166  m_x += m_xInc;
167  m_y += m_yInc;
168  m_error += m_dx;
169  m_error -= m_dy;
170  --m_n;
171  }
172  --m_n;
173  }
174 
175  bool isValid() const { return m_n > 0; }
176 
177  private:
178  float m_dx, m_dy;
179  int m_x, m_y;
180  int m_xInc, m_yInc;
181  float m_error;
182  int m_n;
183 };
184 
186 
187 QgsSnapIndex::GridRow::~GridRow()
188 {
189  Q_FOREACH ( const QgsSnapIndex::Cell &cell, mCells )
190  {
191  qDeleteAll( cell );
192  }
193 }
194 
195 QgsSnapIndex::Cell &QgsSnapIndex::GridRow::getCreateCell( int col )
196 {
197  if ( col < mColStartIdx )
198  {
199  for ( int i = col; i < mColStartIdx; ++i )
200  {
201  mCells.prepend( Cell() );
202  }
203  mColStartIdx = col;
204  return mCells.front();
205  }
206  else if ( col >= mColStartIdx + mCells.size() )
207  {
208  for ( int i = mColStartIdx + mCells.size(); i <= col; ++i )
209  {
210  mCells.append( Cell() );
211  }
212  return mCells.back();
213  }
214  else
215  {
216  return mCells[col - mColStartIdx];
217  }
218 }
219 
220 const QgsSnapIndex::Cell *QgsSnapIndex::GridRow::getCell( int col ) const
221 {
222  if ( col < mColStartIdx || col >= mColStartIdx + mCells.size() )
223  {
224  return nullptr;
225  }
226  else
227  {
228  return &mCells[col - mColStartIdx];
229  }
230 }
231 
232 QList<QgsSnapIndex::SnapItem *> QgsSnapIndex::GridRow::getSnapItems( int colStart, int colEnd ) const
233 {
234  colStart = std::max( colStart, mColStartIdx );
235  colEnd = std::min( colEnd, mColStartIdx + mCells.size() - 1 );
236 
237  QList<SnapItem *> items;
238 
239  for ( int col = colStart; col <= colEnd; ++col )
240  {
241  items.append( mCells[col - mColStartIdx] );
242  }
243  return items;
244 }
245 
247 
248 QgsSnapIndex::QgsSnapIndex( const QgsPoint &origin, double cellSize )
249  : mOrigin( origin )
250  , mCellSize( cellSize )
251  , mRowsStartIdx( 0 )
252 {
253 }
254 
255 QgsSnapIndex::~QgsSnapIndex()
256 {
257  qDeleteAll( mCoordIdxs );
258 }
259 
260 
261 const QgsSnapIndex::Cell *QgsSnapIndex::getCell( int col, int row ) const
262 {
263  if ( row < mRowsStartIdx || row >= mRowsStartIdx + mGridRows.size() )
264  {
265  return nullptr;
266  }
267  else
268  {
269  return mGridRows[row - mRowsStartIdx].getCell( col );
270  }
271 }
272 
273 QgsSnapIndex::Cell &QgsSnapIndex::getCreateCell( int col, int row )
274 {
275  if ( row < mRowsStartIdx )
276  {
277  for ( int i = row; i < mRowsStartIdx; ++i )
278  {
279  mGridRows.prepend( GridRow() );
280  }
281  mRowsStartIdx = row;
282  return mGridRows.front().getCreateCell( col );
283  }
284  else if ( row >= mRowsStartIdx + mGridRows.size() )
285  {
286  for ( int i = mRowsStartIdx + mGridRows.size(); i <= row; ++i )
287  {
288  mGridRows.append( GridRow() );
289  }
290  return mGridRows.back().getCreateCell( col );
291  }
292  else
293  {
294  return mGridRows[row - mRowsStartIdx].getCreateCell( col );
295  }
296 }
297 
298 void QgsSnapIndex::addPoint( const CoordIdx *idx, bool isEndPoint )
299 {
300  QgsPoint p = idx->point();
301  int col = std::floor( ( p.x() - mOrigin.x() ) / mCellSize );
302  int row = std::floor( ( p.y() - mOrigin.y() ) / mCellSize );
303  getCreateCell( col, row ).append( new PointSnapItem( idx, isEndPoint ) );
304 }
305 
306 void QgsSnapIndex::addSegment( const CoordIdx *idxFrom, const CoordIdx *idxTo )
307 {
308  QgsPoint pFrom = idxFrom->point();
309  QgsPoint pTo = idxTo->point();
310  // Raytrace along the grid, get touched cells
311  float x0 = ( pFrom.x() - mOrigin.x() ) / mCellSize;
312  float y0 = ( pFrom.y() - mOrigin.y() ) / mCellSize;
313  float x1 = ( pTo.x() - mOrigin.x() ) / mCellSize;
314  float y1 = ( pTo.y() - mOrigin.y() ) / mCellSize;
315 
316  Raytracer rt( x0, y0, x1, y1 );
317  for ( ; rt.isValid(); rt.next() )
318  {
319  getCreateCell( rt.curCol(), rt.curRow() ).append( new SegmentSnapItem( idxFrom, idxTo ) );
320  }
321 }
322 
323 void QgsSnapIndex::addGeometry( const QgsAbstractGeometry *geom )
324 {
325  for ( int iPart = 0, nParts = geom->partCount(); iPart < nParts; ++iPart )
326  {
327  for ( int iRing = 0, nRings = geom->ringCount( iPart ); iRing < nRings; ++iRing )
328  {
329  int nVerts = geom->vertexCount( iPart, iRing );
330 
331  if ( qgsgeometry_cast< const QgsSurface * >( geom ) )
332  nVerts--;
333  else if ( const QgsCurve *curve = qgsgeometry_cast< const QgsCurve * >( geom ) )
334  {
335  if ( curve->isClosed() )
336  nVerts--;
337  }
338 
339  for ( int iVert = 0; iVert < nVerts; ++iVert )
340  {
341  CoordIdx *idx = new CoordIdx( geom, QgsVertexId( iPart, iRing, iVert ) );
342  CoordIdx *idx1 = new CoordIdx( geom, QgsVertexId( iPart, iRing, iVert + 1 ) );
343  mCoordIdxs.append( idx );
344  mCoordIdxs.append( idx1 );
345  addPoint( idx, iVert == 0 || iVert == nVerts - 1 );
346  if ( iVert < nVerts - 1 )
347  addSegment( idx, idx1 );
348  }
349  }
350  }
351 }
352 
353 
354 QgsPoint QgsSnapIndex::getClosestSnapToPoint( const QgsPoint &p, const QgsPoint &q )
355 {
356  // Look for intersections on segment from the target point to the point opposite to the point reference point
357  // p2 = p1 + 2 * (q - p1)
358  QgsPoint p2( 2 * q.x() - p.x(), 2 * q.y() - p.y() );
359 
360  // Raytrace along the grid, get touched cells
361  float x0 = ( p.x() - mOrigin.x() ) / mCellSize;
362  float y0 = ( p.y() - mOrigin.y() ) / mCellSize;
363  float x1 = ( p2.x() - mOrigin.x() ) / mCellSize;
364  float y1 = ( p2.y() - mOrigin.y() ) / mCellSize;
365 
366  Raytracer rt( x0, y0, x1, y1 );
367  double dMin = std::numeric_limits<double>::max();
368  QgsPoint pMin = p;
369  for ( ; rt.isValid(); rt.next() )
370  {
371  const Cell *cell = getCell( rt.curCol(), rt.curRow() );
372  if ( !cell )
373  {
374  continue;
375  }
376  Q_FOREACH ( const SnapItem *item, *cell )
377  {
378  if ( item->type == SnapSegment )
379  {
380  QgsPoint inter;
381  if ( static_cast<const SegmentSnapItem *>( item )->getIntersection( p, p2, inter ) )
382  {
383  double dist = QgsGeometryUtils::sqrDistance2D( q, inter );
384  if ( dist < dMin )
385  {
386  dMin = dist;
387  pMin = inter;
388  }
389  }
390  }
391  }
392  }
393 
394  return pMin;
395 }
396 
397 QgsSnapIndex::SnapItem *QgsSnapIndex::getSnapItem( const QgsPoint &pos, double tol, QgsSnapIndex::PointSnapItem **pSnapPoint, QgsSnapIndex::SegmentSnapItem **pSnapSegment, bool endPointOnly ) const
398 {
399  int colStart = std::floor( ( pos.x() - tol - mOrigin.x() ) / mCellSize );
400  int rowStart = std::floor( ( pos.y() - tol - mOrigin.y() ) / mCellSize );
401  int colEnd = std::floor( ( pos.x() + tol - mOrigin.x() ) / mCellSize );
402  int rowEnd = std::floor( ( pos.y() + tol - mOrigin.y() ) / mCellSize );
403 
404  rowStart = std::max( rowStart, mRowsStartIdx );
405  rowEnd = std::min( rowEnd, mRowsStartIdx + mGridRows.size() - 1 );
406 
407  QList<SnapItem *> items;
408  for ( int row = rowStart; row <= rowEnd; ++row )
409  {
410  items.append( mGridRows[row - mRowsStartIdx].getSnapItems( colStart, colEnd ) );
411  }
412 
413  double minDistSegment = std::numeric_limits<double>::max();
414  double minDistPoint = std::numeric_limits<double>::max();
415  QgsSnapIndex::SegmentSnapItem *snapSegment = nullptr;
416  QgsSnapIndex::PointSnapItem *snapPoint = nullptr;
417 
418  Q_FOREACH ( QgsSnapIndex::SnapItem *item, items )
419  {
420  if ( ( ! endPointOnly && item->type == SnapPoint ) || item->type == SnapEndPoint )
421  {
422  double dist = QgsGeometryUtils::sqrDistance2D( item->getSnapPoint( pos ), pos );
423  if ( dist < minDistPoint )
424  {
425  minDistPoint = dist;
426  snapPoint = static_cast<PointSnapItem *>( item );
427  }
428  }
429  else if ( item->type == SnapSegment && !endPointOnly )
430  {
431  QgsPoint pProj;
432  if ( !static_cast<SegmentSnapItem *>( item )->getProjection( pos, pProj ) )
433  {
434  continue;
435  }
436  double dist = QgsGeometryUtils::sqrDistance2D( pProj, pos );
437  if ( dist < minDistSegment )
438  {
439  minDistSegment = dist;
440  snapSegment = static_cast<SegmentSnapItem *>( item );
441  }
442  }
443  }
444  snapPoint = minDistPoint < tol * tol ? snapPoint : nullptr;
445  snapSegment = minDistSegment < tol * tol ? snapSegment : nullptr;
446  if ( pSnapPoint ) *pSnapPoint = snapPoint;
447  if ( pSnapSegment ) *pSnapSegment = snapSegment;
448  return minDistPoint < minDistSegment ? static_cast<QgsSnapIndex::SnapItem *>( snapPoint ) : static_cast<QgsSnapIndex::SnapItem *>( snapSegment );
449 }
450 
452 
453 
454 //
455 // QgsGeometrySnapper
456 //
457 
459  : mReferenceSource( referenceSource )
460 {
461  // Build spatial index
462  mIndex = QgsSpatialIndex( *mReferenceSource );
463 }
464 
465 QgsFeatureList QgsGeometrySnapper::snapFeatures( const QgsFeatureList &features, double snapTolerance, SnapMode mode )
466 {
467  QgsFeatureList list = features;
468  QtConcurrent::blockingMap( list, ProcessFeatureWrapper( this, snapTolerance, mode ) );
469  return list;
470 }
471 
472 void QgsGeometrySnapper::processFeature( QgsFeature &feature, double snapTolerance, SnapMode mode )
473 {
474  if ( !feature.geometry().isNull() )
475  feature.setGeometry( snapGeometry( feature.geometry(), snapTolerance, mode ) );
476  emit featureSnapped();
477 }
478 
479 QgsGeometry QgsGeometrySnapper::snapGeometry( const QgsGeometry &geometry, double snapTolerance, SnapMode mode ) const
480 {
481  // Get potential reference features and construct snap index
482  QList<QgsGeometry> refGeometries;
483  mIndexMutex.lock();
484  QgsRectangle searchBounds = geometry.boundingBox();
485  searchBounds.grow( snapTolerance );
486  QgsFeatureIds refFeatureIds = mIndex.intersects( searchBounds ).toSet();
487  mIndexMutex.unlock();
488 
489  QgsFeatureRequest refFeatureRequest = QgsFeatureRequest().setFilterFids( refFeatureIds ).setNoAttributes();
490  mReferenceLayerMutex.lock();
491  QgsFeature refFeature;
492  QgsFeatureIterator refFeatureIt = mReferenceSource->getFeatures( refFeatureRequest );
493 
494  while ( refFeatureIt.nextFeature( refFeature ) )
495  {
496  refGeometries.append( refFeature.geometry() );
497  }
498  mReferenceLayerMutex.unlock();
499 
500  return snapGeometry( geometry, snapTolerance, refGeometries, mode );
501 }
502 
503 QgsGeometry QgsGeometrySnapper::snapGeometry( const QgsGeometry &geometry, double snapTolerance, const QList<QgsGeometry> &referenceGeometries, QgsGeometrySnapper::SnapMode mode )
504 {
506  ( mode == EndPointPreferClosest || mode == EndPointPreferNodes || mode == EndPointToEndPoint ) )
507  return geometry;
508 
509  QgsPoint center = qgsgeometry_cast< const QgsPoint * >( geometry.constGet() ) ? *static_cast< const QgsPoint * >( geometry.constGet() ) :
510  QgsPoint( geometry.constGet()->boundingBox().center() );
511 
512  QgsSnapIndex refSnapIndex( center, 10 * snapTolerance );
513  Q_FOREACH ( const QgsGeometry &geom, referenceGeometries )
514  {
515  refSnapIndex.addGeometry( geom.constGet() );
516  }
517 
518  // Snap geometries
519  QgsAbstractGeometry *subjGeom = geometry.constGet()->clone();
520  QList < QList< QList<PointFlag> > > subjPointFlags;
521 
522  // Pass 1: snap vertices of subject geometry to reference vertices
523  for ( int iPart = 0, nParts = subjGeom->partCount(); iPart < nParts; ++iPart )
524  {
525  subjPointFlags.append( QList< QList<PointFlag> >() );
526 
527  for ( int iRing = 0, nRings = subjGeom->ringCount( iPart ); iRing < nRings; ++iRing )
528  {
529  subjPointFlags[iPart].append( QList<PointFlag>() );
530 
531  for ( int iVert = 0, nVerts = polyLineSize( subjGeom, iPart, iRing ); iVert < nVerts; ++iVert )
532  {
533  if ( ( mode == EndPointPreferClosest || mode == EndPointPreferNodes || mode == EndPointToEndPoint ) &&
534  QgsWkbTypes::geometryType( subjGeom->wkbType() ) == QgsWkbTypes::LineGeometry && ( iVert > 0 && iVert < nVerts - 1 ) )
535  {
536  //endpoint mode and not at an endpoint, skip
537  subjPointFlags[iPart][iRing].append( Unsnapped );
538  continue;
539  }
540 
541  QgsSnapIndex::PointSnapItem *snapPoint = nullptr;
542  QgsSnapIndex::SegmentSnapItem *snapSegment = nullptr;
543  QgsVertexId vidx( iPart, iRing, iVert );
544  QgsPoint p = subjGeom->vertexAt( vidx );
545  if ( !refSnapIndex.getSnapItem( p, snapTolerance, &snapPoint, &snapSegment, mode == EndPointToEndPoint ) )
546  {
547  subjPointFlags[iPart][iRing].append( Unsnapped );
548  }
549  else
550  {
551  switch ( mode )
552  {
553  case PreferNodes:
555  case EndPointPreferNodes:
556  case EndPointToEndPoint:
557  {
558  // Prefer snapping to point
559  if ( snapPoint )
560  {
561  subjGeom->moveVertex( vidx, snapPoint->getSnapPoint( p ) );
562  subjPointFlags[iPart][iRing].append( SnappedToRefNode );
563  }
564  else if ( snapSegment )
565  {
566  subjGeom->moveVertex( vidx, snapSegment->getSnapPoint( p ) );
567  subjPointFlags[iPart][iRing].append( SnappedToRefSegment );
568  }
569  break;
570  }
571 
572  case PreferClosest:
575  {
576  QgsPoint nodeSnap, segmentSnap;
577  double distanceNode = std::numeric_limits<double>::max();
578  double distanceSegment = std::numeric_limits<double>::max();
579  if ( snapPoint )
580  {
581  nodeSnap = snapPoint->getSnapPoint( p );
582  distanceNode = nodeSnap.distanceSquared( p );
583  }
584  if ( snapSegment )
585  {
586  segmentSnap = snapSegment->getSnapPoint( p );
587  distanceSegment = segmentSnap.distanceSquared( p );
588  }
589  if ( snapPoint && distanceNode < distanceSegment )
590  {
591  subjGeom->moveVertex( vidx, nodeSnap );
592  subjPointFlags[iPart][iRing].append( SnappedToRefNode );
593  }
594  else if ( snapSegment )
595  {
596  subjGeom->moveVertex( vidx, segmentSnap );
597  subjPointFlags[iPart][iRing].append( SnappedToRefSegment );
598  }
599  break;
600  }
601  }
602  }
603  }
604  }
605  }
606 
607  //nothing more to do for points
608  if ( qgsgeometry_cast< const QgsPoint * >( subjGeom ) )
609  return QgsGeometry( subjGeom );
610  //or for end point snapping
612  {
613  QgsGeometry result( subjGeom );
614  result.removeDuplicateNodes();
615  return result;
616  }
617 
618  // SnapIndex for subject feature
619  std::unique_ptr< QgsSnapIndex > subjSnapIndex( new QgsSnapIndex( center, 10 * snapTolerance ) );
620  subjSnapIndex->addGeometry( subjGeom );
621 
622  std::unique_ptr< QgsAbstractGeometry > origSubjGeom( subjGeom->clone() );
623  std::unique_ptr< QgsSnapIndex > origSubjSnapIndex( new QgsSnapIndex( center, 10 * snapTolerance ) );
624  origSubjSnapIndex->addGeometry( origSubjGeom.get() );
625 
626  // Pass 2: add missing vertices to subject geometry
627  Q_FOREACH ( const QgsGeometry &refGeom, referenceGeometries )
628  {
629  for ( int iPart = 0, nParts = refGeom.constGet()->partCount(); iPart < nParts; ++iPart )
630  {
631  for ( int iRing = 0, nRings = refGeom.constGet()->ringCount( iPart ); iRing < nRings; ++iRing )
632  {
633  for ( int iVert = 0, nVerts = polyLineSize( refGeom.constGet(), iPart, iRing ); iVert < nVerts; ++iVert )
634  {
635 
636  QgsSnapIndex::PointSnapItem *snapPoint = nullptr;
637  QgsSnapIndex::SegmentSnapItem *snapSegment = nullptr;
638  QgsPoint point = refGeom.constGet()->vertexAt( QgsVertexId( iPart, iRing, iVert ) );
639  if ( subjSnapIndex->getSnapItem( point, snapTolerance, &snapPoint, &snapSegment ) )
640  {
641  // Snap to segment, unless a subject point was already snapped to the reference point
642  if ( snapPoint && QgsGeometryUtils::sqrDistance2D( snapPoint->getSnapPoint( point ), point ) < 1E-16 )
643  {
644  continue;
645  }
646  else if ( snapSegment )
647  {
648  // Look if there is a closer reference segment, if so, ignore this point
649  QgsPoint pProj = snapSegment->getSnapPoint( point );
650  QgsPoint closest = refSnapIndex.getClosestSnapToPoint( point, pProj );
651  if ( QgsGeometryUtils::sqrDistance2D( pProj, point ) > QgsGeometryUtils::sqrDistance2D( pProj, closest ) )
652  {
653  continue;
654  }
655 
656  // If we are too far away from the original geometry, do nothing
657  if ( !origSubjSnapIndex->getSnapItem( point, snapTolerance ) )
658  {
659  continue;
660  }
661 
662  const QgsSnapIndex::CoordIdx *idx = snapSegment->idxFrom;
663  subjGeom->insertVertex( QgsVertexId( idx->vidx.part, idx->vidx.ring, idx->vidx.vertex + 1 ), point );
664  subjPointFlags[idx->vidx.part][idx->vidx.ring].insert( idx->vidx.vertex + 1, SnappedToRefNode );
665  subjSnapIndex.reset( new QgsSnapIndex( center, 10 * snapTolerance ) );
666  subjSnapIndex->addGeometry( subjGeom );
667  }
668  }
669  }
670  }
671  }
672  }
673  subjSnapIndex.reset();
674  origSubjSnapIndex.reset();
675  origSubjGeom.reset();
676 
677  // Pass 3: remove superfluous vertices: all vertices which are snapped to a segment and not preceded or succeeded by an unsnapped vertex
678  for ( int iPart = 0, nParts = subjGeom->partCount(); iPart < nParts; ++iPart )
679  {
680  for ( int iRing = 0, nRings = subjGeom->ringCount( iPart ); iRing < nRings; ++iRing )
681  {
682  bool ringIsClosed = subjGeom->vertexAt( QgsVertexId( iPart, iRing, 0 ) ) == subjGeom->vertexAt( QgsVertexId( iPart, iRing, subjGeom->vertexCount( iPart, iRing ) - 1 ) );
683  for ( int iVert = 0, nVerts = polyLineSize( subjGeom, iPart, iRing ); iVert < nVerts; ++iVert )
684  {
685  int iPrev = ( iVert - 1 + nVerts ) % nVerts;
686  int iNext = ( iVert + 1 ) % nVerts;
687  QgsPoint pMid = subjGeom->vertexAt( QgsVertexId( iPart, iRing, iVert ) );
688  QgsPoint pPrev = subjGeom->vertexAt( QgsVertexId( iPart, iRing, iPrev ) );
689  QgsPoint pNext = subjGeom->vertexAt( QgsVertexId( iPart, iRing, iNext ) );
690 
691  if ( subjPointFlags[iPart][iRing][iVert] == SnappedToRefSegment &&
692  subjPointFlags[iPart][iRing][iPrev] != Unsnapped &&
693  subjPointFlags[iPart][iRing][iNext] != Unsnapped &&
694  QgsGeometryUtils::sqrDistance2D( QgsGeometryUtils::projectPointOnSegment( pMid, pPrev, pNext ), pMid ) < 1E-12 )
695  {
696  if ( ( ringIsClosed && nVerts > 3 ) || ( !ringIsClosed && nVerts > 2 ) )
697  {
698  subjGeom->deleteVertex( QgsVertexId( iPart, iRing, iVert ) );
699  subjPointFlags[iPart][iRing].removeAt( iVert );
700  iVert -= 1;
701  nVerts -= 1;
702  }
703  else
704  {
705  // Don't delete vertices if this would result in a degenerate geometry
706  break;
707  }
708  }
709  }
710  }
711  }
712 
713  QgsGeometry result( subjGeom );
714  result.removeDuplicateNodes();
715  return result;
716 }
717 
718 int QgsGeometrySnapper::polyLineSize( const QgsAbstractGeometry *geom, int iPart, int iRing )
719 {
720  int nVerts = geom->vertexCount( iPart, iRing );
721 
722  if ( qgsgeometry_cast< const QgsSurface * >( geom ) || qgsgeometry_cast< const QgsMultiSurface * >( geom ) )
723  {
724  QgsPoint front = geom->vertexAt( QgsVertexId( iPart, iRing, 0 ) );
725  QgsPoint back = geom->vertexAt( QgsVertexId( iPart, iRing, nVerts - 1 ) );
726  if ( front == back )
727  return nVerts - 1;
728  }
729 
730  return nVerts;
731 }
732 
733 
734 
735 
736 
737 //
738 // QgsInternalGeometrySnapper
739 //
740 
742  : mSnapTolerance( snapTolerance )
743  , mMode( mode )
744 {}
745 
747 {
748  if ( !feature.hasGeometry() )
749  return QgsGeometry();
750 
751  QgsFeature feat = feature;
752  QgsGeometry geometry = feat.geometry();
753  if ( !mFirstFeature )
754  {
755  // snap against processed geometries
756  // Get potential reference features and construct snap index
757  QgsRectangle searchBounds = geometry.boundingBox();
758  searchBounds.grow( mSnapTolerance );
759  QgsFeatureIds refFeatureIds = mProcessedIndex.intersects( searchBounds ).toSet();
760  if ( !refFeatureIds.isEmpty() )
761  {
762  QList< QgsGeometry > refGeometries;
763  Q_FOREACH ( QgsFeatureId id, refFeatureIds )
764  {
765  refGeometries << mProcessedGeometries.value( id );
766  }
767 
768  geometry = QgsGeometrySnapper::snapGeometry( geometry, mSnapTolerance, refGeometries, mMode );
769  }
770  }
771  mProcessedGeometries.insert( feat.id(), geometry );
772  mProcessedIndex.addFeature( feat );
773  mFirstFeature = false;
774  return geometry;
775 }
QgsFeatureId id
Definition: qgsfeature.h:64
Wrapper for iterator of features from vector data provider or vector layer.
A rectangle specified with double values.
Definition: qgsrectangle.h:41
double y
Definition: qgspoint.h:42
QSet< QgsFeatureId > QgsFeatureIds
Definition: qgsfeatureid.h:34
virtual bool deleteVertex(QgsVertexId position)=0
Deletes a vertex within the geometry.
bool removeDuplicateNodes(double epsilon=4 *std::numeric_limits< double >::epsilon(), bool useZValues=false)
Removes duplicate nodes from the geometry, wherever removing the nodes does not result in a degenerat...
virtual bool insertVertex(QgsVertexId position, const QgsPoint &vertex)=0
Inserts a vertex into the geometry.
QList< QgsFeature > QgsFeatureList
Definition: qgsfeature.h:571
QgsWkbTypes::Type wkbType() const
Returns type of the geometry as a WKB type (point / linestring / polygon etc.)
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
Only snap the start/end points of lines to other start/end points of lines.
QgsGeometrySnapper(QgsFeatureSource *referenceSource)
Constructor for QgsGeometrySnapper.
qint64 QgsFeatureId
Definition: qgsfeatureid.h:25
void featureSnapped()
Emitted each time a feature has been processed when calling snapFeatures()
A geometry is the spatial representation of a feature.
Definition: qgsgeometry.h:106
The feature class encapsulates a single feature including its id, geometry and a list of field/values...
Definition: qgsfeature.h:55
virtual QgsRectangle boundingBox() const =0
Returns the minimal bounding box for the geometry.
bool hasGeometry() const
Returns true if the feature has an associated geometry.
Definition: qgsfeature.cpp:197
virtual QgsAbstractGeometry * clone() const =0
Clones the geometry by performing a deep copy.
QgsFeatureList snapFeatures(const QgsFeatureList &features, double snapTolerance, SnapMode mode=PreferNodes)
Snaps a set of features to the reference layer and returns the result.
Utility class for identifying a unique vertex within a geometry.
void grow(double delta)
Grows the rectangle in place by the specified amount.
Definition: qgsrectangle.h:275
Only snap start/end points of lines (point features will also be snapped, polygon features will not b...
QgsFeatureRequest & setNoAttributes()
Set that no attributes will be fetched.
static GeometryType geometryType(Type type)
Returns the geometry type for a WKB type, e.g., both MultiPolygon and CurvePolygon would have a Polyg...
Definition: qgswkbtypes.h:665
This class wraps a request for features to a vector layer (or directly its vector data provider)...
T qgsgeometry_cast(const QgsAbstractGeometry *geom)
virtual int ringCount(int part=0) const =0
Returns the number of rings of which this geometry is built.
Abstract base class for curved geometry type.
Definition: qgscurve.h:35
Abstract base class for all geometries.
QgsWkbTypes::Type wkbType() const
Returns the WKB type of the geometry.
Point geometry type, with support for z-dimension and m-values.
Definition: qgspoint.h:37
double length() const
Returns the length of the vector.
Definition: qgsvector.cpp:71
A class to represent a vector.
Definition: qgsvector.h:29
const QgsAbstractGeometry * constGet() const
Returns a non-modifiable (const) reference to the underlying abstract geometry primitive.
QgsGeometry snapGeometry(const QgsGeometry &geometry, double snapTolerance, SnapMode mode=PreferNodes) const
Snaps a geometry to the reference layer and returns the result.
Prefer to snap to nodes, even when a segment may be closer than a node. New nodes will be inserted to...
static double sqrDistance2D(const QgsPoint &pt1, const QgsPoint &pt2)
Returns the squared 2D distance between two points.
A spatial index for QgsFeature objects.
QgsFeatureRequest & setFilterFids(const QgsFeatureIds &fids)
Sets feature IDs that should be fetched.
Snap to closest point, regardless of it is a node or a segment. New nodes will be inserted to make ge...
An interface for objects which provide features via a getFeatures method.
SnapMode
Snapping modes.
virtual bool moveVertex(QgsVertexId position, const QgsPoint &newPos)=0
Moves a vertex within the geometry.
QgsGeometry snapFeature(const QgsFeature &feature)
Snaps a single feature&#39;s geometry against all feature geometries already processed by calls to snapFe...
virtual int vertexCount(int part=0, int ring=0) const =0
Returns the number of vertices of which this geometry is built.
QgsRectangle boundingBox() const
Returns the bounding box of the geometry.
void setGeometry(const QgsGeometry &geometry)
Set the feature&#39;s geometry.
Definition: qgsfeature.cpp:137
static QgsPoint projectPointOnSegment(const QgsPoint &p, const QgsPoint &s1, const QgsPoint &s2)
Project the point on a segment.
QgsInternalGeometrySnapper(double snapTolerance, QgsGeometrySnapper::SnapMode mode=QgsGeometrySnapper::PreferNodes)
Constructor for QgsInternalGeometrySnapper.
Prefer to snap to nodes, even when a segment may be closer than a node. No new nodes will be inserted...
QgsGeometry geometry
Definition: qgsfeature.h:67
QgsPointXY center() const
Returns the center point of the rectangle.
Definition: qgsrectangle.h:230
Only snap start/end points of lines (point features will also be snapped, polygon features will not b...
bool addFeature(QgsFeature &feature, QgsFeatureSink::Flags flags=nullptr) override
Adds a feature to the index.
bool nextFeature(QgsFeature &f)
double distanceSquared(double x, double y) const
Returns the squared distance between this point a specified x, y coordinate.
Definition: qgspoint.h:299
virtual QgsPoint vertexAt(QgsVertexId id) const =0
Returns the point corresponding to a specified vertex id.
Snap to closest point, regardless of it is a node or a segment. No new nodes will be inserted...
QList< QgsFeatureId > intersects(const QgsRectangle &rectangle) const
Returns a list of features with a bounding box which intersects the specified rectangle.
virtual int partCount() const =0
Returns count of parts contained in the geometry.
virtual QgsFeatureIterator getFeatures(const QgsFeatureRequest &request=QgsFeatureRequest()) const =0
Returns an iterator for the features in the source.
double x
Definition: qgspoint.h:41