27#include <QtConcurrentMap>
32QgsSnapIndex::PointSnapItem::PointSnapItem(
const QgsSnapIndex::CoordIdx *_idx,
bool isEndPoint )
33 : SnapItem( isEndPoint ? QgsSnapIndex::SnapEndPoint : QgsSnapIndex::SnapPoint )
42QgsSnapIndex::SegmentSnapItem::SegmentSnapItem(
const QgsSnapIndex::CoordIdx *_idxFrom,
const QgsSnapIndex::CoordIdx *_idxTo )
43 : SnapItem( QgsSnapIndex::SnapSegment )
48QgsPoint QgsSnapIndex::SegmentSnapItem::getSnapPoint(
const QgsPoint &p )
const
55 const QgsPoint &q1 = idxFrom->point(), & q2 = idxTo->point();
58 const double vl = v.length();
59 const double wl = w.length();
68 const double d = v.y() * w.x() - v.x() * w.y();
73 const double dx = q1.
x() - p1.
x();
74 const double dy = q1.
y() - p1.
y();
75 const double k = ( dy * w.x() - dx * w.y() ) / d;
77 inter =
QgsPoint( p1.
x() + v.x() * k, p1.
y() + v.y() * k );
79 const double lambdav =
QgsVector( inter.
x() - p1.
x(), inter.
y() - p1.
y() ) * v;
80 if ( lambdav < 0. + 1E-8 || lambdav > vl - 1E-8 )
83 const double lambdaw =
QgsVector( inter.
x() - q1.
x(), inter.
y() - q1.
y() ) * w;
84 return !( lambdaw < 0. + 1E-8 || lambdaw >= wl - 1E-8 );
87bool QgsSnapIndex::SegmentSnapItem::getProjection(
const QgsPoint &p,
QgsPoint &pProj )
const
89 const QgsPoint &s1 = idxFrom->point();
91 const double nx = s2.
y() - s1.
y();
92 const double ny = -( s2.
x() - s1.
x() );
93 const double t = ( p.
x() * ny - p.
y() * nx - s1.
x() * ny + s1.
y() * nx ) / ( ( s2.
x() - s1.
x() ) * ny - ( s2.
y() - s1.
y() ) * nx );
94 if ( t < 0. || t > 1. )
98 pProj =
QgsPoint( s1.
x() + ( s2.
x() - s1.
x() ) * t, s1.
y() + ( s2.
y() - s1.
y() ) * t );
102bool QgsSnapIndex::SegmentSnapItem::withinSquaredDistance(
const QgsPoint &p,
const double squaredDistance )
104 double minDistX, minDistY;
105 return QgsGeometryUtilsBase::sqrDistToLine( p.
x(), p.
y(), idxFrom->point().x(), idxFrom->point().y(), idxTo->point().x(), idxTo->point().y(), minDistX, minDistY, 4 * std::numeric_limits<double>::epsilon() ) <= squaredDistance;
110QgsSnapIndex::QgsSnapIndex()
115QgsSnapIndex::~QgsSnapIndex()
117 qDeleteAll( mCoordIdxs );
118 qDeleteAll( mSnapItems );
123void QgsSnapIndex::addPoint(
const CoordIdx *idx,
bool isEndPoint )
128 geos::unique_ptr point( GEOSGeom_createPointFromXY_r( geosctxt, p.
x(), p.
y() ) );
130 PointSnapItem *item =
new PointSnapItem( idx, isEndPoint );
131 GEOSSTRtree_insert_r( geosctxt, mSTRTree, point.get(), item );
135void QgsSnapIndex::addSegment(
const CoordIdx *idxFrom,
const CoordIdx *idxTo )
137 const QgsPoint pointFrom = idxFrom->point();
138 const QgsPoint pointTo = idxTo->point();
142 GEOSCoordSequence *coord = GEOSCoordSeq_create_r( geosctxt, 2, 2 );
143 GEOSCoordSeq_setXY_r( geosctxt, coord, 0, pointFrom.
x(), pointFrom.
y() );
144 GEOSCoordSeq_setXY_r( geosctxt, coord, 1, pointTo.
x(), pointTo.
y() );
145 geos::unique_ptr
segment( GEOSGeom_createLineString_r( geosctxt, coord ) );
147 SegmentSnapItem *item =
new SegmentSnapItem( idxFrom, idxTo );
148 GEOSSTRtree_insert_r( geosctxt, mSTRTree,
segment.get(), item );
154 for (
int iPart = 0, nParts = geom->
partCount(); iPart < nParts; ++iPart )
156 for (
int iRing = 0, nRings = geom->
ringCount( iPart ); iRing < nRings; ++iRing )
160 if ( qgsgeometry_cast< const QgsSurface * >( geom ) )
162 else if (
const QgsCurve *curve = qgsgeometry_cast< const QgsCurve * >( geom ) )
164 if ( curve->isClosed() )
168 for (
int iVert = 0; iVert < nVerts; ++iVert )
170 CoordIdx *idx =
new CoordIdx( geom,
QgsVertexId( iPart, iRing, iVert ) );
171 CoordIdx *idx1 =
new CoordIdx( geom,
QgsVertexId( iPart, iRing, iVert + 1 ) );
172 mCoordIdxs.append( idx );
173 mCoordIdxs.append( idx1 );
174 addPoint( idx, iVert == 0 || iVert == nVerts - 1 );
175 if ( iVert < nVerts - 1 )
176 addSegment( idx, idx1 );
184 QList< QgsSnapIndex::SnapItem * > *
list;
189 reinterpret_cast<_GEOSQueryCallbackData *
>( userdata )->list->append(
static_cast<QgsSnapIndex::SnapItem *
>( item ) );
198 const QgsPoint endPoint( 2 * midPoint.
x() - startPoint.
x(), 2 * midPoint.
y() - startPoint.
y() );
201 double minDistance = std::numeric_limits<double>::max();
203 GEOSCoordSequence *coord = GEOSCoordSeq_create_r( geosctxt, 2, 2 );
204 GEOSCoordSeq_setXY_r( geosctxt, coord, 0, startPoint.
x(), startPoint.
y() );
205 GEOSCoordSeq_setXY_r( geosctxt, coord, 1, endPoint.x(), endPoint.y() );
206 geos::unique_ptr searchDiagonal( GEOSGeom_createLineString_r( geosctxt, coord ) );
208 QList<SnapItem *> items;
210 callbackData.
list = &items;
211 GEOSSTRtree_query_r( geosctxt, mSTRTree, searchDiagonal.get(),
_GEOSQueryCallback, &callbackData );
212 for (
const SnapItem *item : std::as_const( items ) )
214 if ( item->type == SnapSegment )
217 if (
static_cast<const SegmentSnapItem *
>( item )->getIntersection( startPoint, endPoint, inter ) )
220 if ( dist < minDistance )
232QgsSnapIndex::SnapItem *QgsSnapIndex::getSnapItem(
const QgsPoint &pos,
const double tolerance, QgsSnapIndex::PointSnapItem **pSnapPoint, QgsSnapIndex::SegmentSnapItem **pSnapSegment,
bool endPointOnly )
const
236 GEOSCoordSequence *coord = GEOSCoordSeq_create_r( geosctxt, 2, 2 );
237 GEOSCoordSeq_setXY_r( geosctxt, coord, 0, pos.
x() - tolerance, pos.
y() - tolerance );
238 GEOSCoordSeq_setXY_r( geosctxt, coord, 1, pos.
x() + tolerance, pos.
y() + tolerance );
240 geos::unique_ptr searchDiagonal( GEOSGeom_createLineString_r( geosctxt, coord ) );
242 QList<SnapItem *> items;
244 callbackData.
list = &items;
245 GEOSSTRtree_query_r( geosctxt, mSTRTree, searchDiagonal.get(),
_GEOSQueryCallback, &callbackData );
247 double minDistSegment = std::numeric_limits<double>::max();
248 double minDistPoint = std::numeric_limits<double>::max();
249 QgsSnapIndex::SegmentSnapItem *snapSegment =
nullptr;
250 QgsSnapIndex::PointSnapItem *snapPoint =
nullptr;
252 const double squaredTolerance = tolerance * tolerance;
253 const auto constItems = items;
254 for ( QgsSnapIndex::SnapItem *item : constItems )
256 if ( ( ! endPointOnly && item->type == SnapPoint ) || item->type == SnapEndPoint )
259 if ( dist < minDistPoint )
262 snapPoint =
static_cast<PointSnapItem *
>( item );
265 else if ( item->type == SnapSegment && !endPointOnly )
267 if ( !
static_cast<SegmentSnapItem *
>( item )->withinSquaredDistance( pos, squaredTolerance ) )
271 if ( !
static_cast<SegmentSnapItem *
>( item )->getProjection( pos, pProj ) )
275 if ( dist < minDistSegment )
277 minDistSegment = dist;
278 snapSegment =
static_cast<SegmentSnapItem *
>( item );
282 snapPoint = minDistPoint < squaredTolerance ? snapPoint :
nullptr;
283 snapSegment = minDistSegment < squaredTolerance ? snapSegment :
nullptr;
284 if ( pSnapPoint ) *pSnapPoint = snapPoint;
285 if ( pSnapSegment ) *pSnapSegment = snapSegment;
286 return minDistPoint < minDistSegment ? static_cast<QgsSnapIndex::SnapItem *>( snapPoint ) : static_cast<QgsSnapIndex::SnapItem *>( snapSegment );
297 : mReferenceSource( referenceSource )
306 QtConcurrent::blockingMap( list, ProcessFeatureWrapper(
this, snapTolerance, mode ) );
310void QgsGeometrySnapper::processFeature(
QgsFeature &feature,
double snapTolerance, SnapMode mode )
320 QList<QgsGeometry> refGeometries;
323 searchBounds.
grow( snapTolerance );
325 mIndexMutex.unlock();
327 if ( refFeatureIds.isEmpty() )
330 refGeometries.reserve( refFeatureIds.size() );
332 const QgsFeatureIds cachedIds = qgis::listToSet( mCachedReferenceGeometries.keys() );
335 if ( cachedIds.contains(
id ) )
337 refGeometries.append( mCachedReferenceGeometries[
id] );
341 missingFeatureIds << id;
345 if ( missingFeatureIds.size() > 0 )
348 mReferenceLayerMutex.lock();
354 refGeometries.append( refFeature.
geometry() );
356 mReferenceLayerMutex.unlock();
359 return snapGeometry( geometry, snapTolerance, refGeometries, mode );
371 QgsSnapIndex refSnapIndex;
372 for (
const QgsGeometry &geom : referenceGeometries )
374 refSnapIndex.addGeometry( geom.constGet() );
379 QList < QList< QList<PointFlag> > > subjPointFlags;
382 for (
int iPart = 0, nParts = subjGeom->
partCount(); iPart < nParts; ++iPart )
384 subjPointFlags.append( QList< QList<PointFlag> >() );
386 for (
int iRing = 0, nRings = subjGeom->
ringCount( iPart ); iRing < nRings; ++iRing )
388 subjPointFlags[iPart].append( QList<PointFlag>() );
390 for (
int iVert = 0, nVerts = polyLineSize( subjGeom, iPart, iRing ); iVert < nVerts; ++iVert )
396 subjPointFlags[iPart][iRing].append( Unsnapped );
400 QgsSnapIndex::PointSnapItem *snapPoint =
nullptr;
401 QgsSnapIndex::SegmentSnapItem *snapSegment =
nullptr;
404 if ( !refSnapIndex.getSnapItem( p, snapTolerance, &snapPoint, &snapSegment, mode ==
EndPointToEndPoint ) )
406 subjPointFlags[iPart][iRing].append( Unsnapped );
420 subjGeom->
moveVertex( vidx, snapPoint->getSnapPoint( p ) );
421 subjPointFlags[iPart][iRing].append( SnappedToRefNode );
423 else if ( snapSegment )
425 subjGeom->
moveVertex( vidx, snapSegment->getSnapPoint( p ) );
426 subjPointFlags[iPart][iRing].append( SnappedToRefSegment );
436 double distanceNode = std::numeric_limits<double>::max();
437 double distanceSegment = std::numeric_limits<double>::max();
440 nodeSnap = snapPoint->getSnapPoint( p );
445 segmentSnap = snapSegment->getSnapPoint( p );
448 if ( snapPoint && distanceNode < distanceSegment )
451 subjPointFlags[iPart][iRing].append( SnappedToRefNode );
453 else if ( snapSegment )
456 subjPointFlags[iPart][iRing].append( SnappedToRefSegment );
467 if ( qgsgeometry_cast< const QgsPoint * >( subjGeom ) )
478 std::unique_ptr< QgsSnapIndex > subjSnapIndex(
new QgsSnapIndex() );
479 subjSnapIndex->addGeometry( subjGeom );
481 std::unique_ptr< QgsAbstractGeometry > origSubjGeom( subjGeom->
clone() );
482 std::unique_ptr< QgsSnapIndex > origSubjSnapIndex(
new QgsSnapIndex() );
483 origSubjSnapIndex->addGeometry( origSubjGeom.get() );
486 for (
const QgsGeometry &refGeom : referenceGeometries )
488 for (
int iPart = 0, nParts = refGeom.constGet()->partCount(); iPart < nParts; ++iPart )
490 for (
int iRing = 0, nRings = refGeom.constGet()->ringCount( iPart ); iRing < nRings; ++iRing )
492 for (
int iVert = 0, nVerts = polyLineSize( refGeom.constGet(), iPart, iRing ); iVert < nVerts; ++iVert )
494 QgsSnapIndex::PointSnapItem *snapPoint =
nullptr;
495 QgsSnapIndex::SegmentSnapItem *snapSegment =
nullptr;
497 if ( subjSnapIndex->getSnapItem( point, snapTolerance, &snapPoint, &snapSegment ) )
502 const QgsPoint snappedPoint = snapPoint->getSnapPoint( point );
510 const QgsPoint pProj = snapSegment->getSnapPoint( point );
511 const QgsPoint closest = refSnapIndex.getClosestSnapToPoint( point, pProj );
518 if ( !origSubjSnapIndex->getSnapItem( point, snapTolerance ) )
523 const QgsSnapIndex::CoordIdx *idx = snapSegment->idxFrom;
525 subjPointFlags[idx->vidx.part][idx->vidx.ring].insert( idx->vidx.vertex + 1, SnappedToRefNode );
526 subjSnapIndex.reset(
new QgsSnapIndex() );
527 subjSnapIndex->addGeometry( subjGeom );
534 subjSnapIndex.reset();
535 origSubjSnapIndex.reset();
536 origSubjGeom.reset();
539 for (
int iPart = 0, nParts = subjGeom->
partCount(); iPart < nParts; ++iPart )
541 for (
int iRing = 0, nRings = subjGeom->
ringCount( iPart ); iRing < nRings; ++iRing )
544 for (
int iVert = 0, nVerts = polyLineSize( subjGeom, iPart, iRing ); iVert < nVerts; ++iVert )
546 const int iPrev = ( iVert - 1 + nVerts ) % nVerts;
547 const int iNext = ( iVert + 1 ) % nVerts;
552 if ( subjPointFlags[iPart][iRing][iVert] == SnappedToRefSegment &&
553 subjPointFlags[iPart][iRing][iPrev] != Unsnapped &&
554 subjPointFlags[iPart][iRing][iNext] != Unsnapped &&
557 if ( ( ringIsClosed && nVerts > 3 ) || ( !ringIsClosed && nVerts > 2 ) )
560 subjPointFlags[iPart][iRing].removeAt( iVert );
579int QgsGeometrySnapper::polyLineSize(
const QgsAbstractGeometry *geom,
int iPart,
int iRing )
581 const int nVerts = geom->
vertexCount( iPart, iRing );
583 if ( qgsgeometry_cast< const QgsSurface * >( geom ) || qgsgeometry_cast< const QgsMultiSurface * >( geom ) )
603 : mSnapTolerance( snapTolerance )
614 if ( !mFirstFeature )
619 searchBounds.
grow( mSnapTolerance );
621 if ( !refFeatureIds.isEmpty() )
623 QList< QgsGeometry > refGeometries;
624 const auto constRefFeatureIds = refFeatureIds;
627 refGeometries << mProcessedGeometries.value(
id );
633 mProcessedGeometries.insert( feat.
id(), geometry );
635 mFirstFeature =
false;
Abstract base class for all geometries.
virtual int ringCount(int part=0) const =0
Returns the number of rings of which this geometry is built.
virtual bool moveVertex(QgsVertexId position, const QgsPoint &newPos)=0
Moves a vertex within the geometry.
virtual int vertexCount(int part=0, int ring=0) const =0
Returns the number of vertices of which this geometry is built.
virtual QgsRectangle boundingBox() const
Returns the minimal bounding box for the geometry.
virtual QgsPoint vertexAt(QgsVertexId id) const =0
Returns the point corresponding to a specified vertex id.
Qgis::WkbType wkbType() const
Returns the WKB type of the geometry.
virtual bool insertVertex(QgsVertexId position, const QgsPoint &vertex)=0
Inserts a vertex into the geometry.
virtual int partCount() const =0
Returns count of parts contained in the geometry.
virtual bool deleteVertex(QgsVertexId position)=0
Deletes a vertex within the geometry.
virtual QgsAbstractGeometry * clone() const =0
Clones the geometry by performing a deep copy.
Abstract base class for curved geometry type.
Wrapper for iterator of features from vector data provider or vector layer.
bool nextFeature(QgsFeature &f)
Fetch next feature and stores in f, returns true on success.
This class wraps a request for features to a vector layer (or directly its vector data provider).
QgsFeatureRequest & setFilterFids(const QgsFeatureIds &fids)
Sets the feature IDs that should be fetched.
QgsFeatureRequest & setNoAttributes()
Set that no attributes will be fetched.
An interface for objects which provide features via a getFeatures method.
virtual QgsFeatureIterator getFeatures(const QgsFeatureRequest &request=QgsFeatureRequest()) const =0
Returns an iterator for the features in the source.
The feature class encapsulates a single feature including its unique ID, geometry and a list of field...
bool hasGeometry() const
Returns true if the feature has an associated geometry.
void setGeometry(const QgsGeometry &geometry)
Set the feature's geometry.
void featureSnapped()
Emitted each time a feature has been processed when calling snapFeatures()
QgsFeatureList snapFeatures(const QgsFeatureList &features, double snapTolerance, SnapMode mode=PreferNodes)
Snaps a set of features to the reference layer and returns the result.
QgsGeometry snapGeometry(const QgsGeometry &geometry, double snapTolerance, SnapMode mode=PreferNodes) const
Snaps a geometry to the reference layer and returns the result.
@ EndPointPreferClosest
Only snap start/end points of lines (point features will also be snapped, polygon features will not b...
@ PreferClosestNoExtraVertices
Snap to closest point, regardless of it is a node or a segment. No new nodes will be inserted.
@ EndPointPreferNodes
Only snap start/end points of lines (point features will also be snapped, polygon features will not b...
@ PreferNodes
Prefer to snap to nodes, even when a segment may be closer than a node. New nodes will be inserted to...
@ PreferClosest
Snap to closest point, regardless of it is a node or a segment. New nodes will be inserted to make ge...
@ EndPointToEndPoint
Only snap the start/end points of lines to other start/end points of lines.
@ PreferNodesNoExtraVertices
Prefer to snap to nodes, even when a segment may be closer than a node. No new nodes will be inserted...
QgsGeometrySnapper(QgsFeatureSource *referenceSource)
Constructor for QgsGeometrySnapper.
static double sqrDistToLine(double ptX, double ptY, double x1, double y1, double x2, double y2, double &minDistX, double &minDistY, double epsilon)
Returns the squared distance between a point and a line.
static QgsPoint projectPointOnSegment(const QgsPoint &p, const QgsPoint &s1, const QgsPoint &s2)
Project the point on a segment.
static Q_DECL_DEPRECATED double sqrDistance2D(double x1, double y1, double x2, double y2)
Returns the squared 2D distance between (x1, y1) and (x2, y2).
A geometry is the spatial representation of a feature.
const QgsAbstractGeometry * constGet() const
Returns a non-modifiable (const) reference to the underlying abstract geometry primitive.
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...
QgsRectangle boundingBox() const
Returns the bounding box of the geometry.
Qgis::WkbType wkbType() const
Returns type of the geometry as a WKB type (point / linestring / polygon etc.)
static GEOSContextHandle_t get()
Returns a thread local instance of a GEOS context, safe for use in the current thread.
QgsInternalGeometrySnapper(double snapTolerance, QgsGeometrySnapper::SnapMode mode=QgsGeometrySnapper::PreferNodes)
Constructor for QgsInternalGeometrySnapper.
QgsGeometry snapFeature(const QgsFeature &feature)
Snaps a single feature's geometry against all feature geometries already processed by calls to snapFe...
Point geometry type, with support for z-dimension and m-values.
QgsPoint vertexAt(QgsVertexId) const override
Returns the point corresponding to a specified vertex id.
double distanceSquared(double x, double y) const
Returns the Cartesian 2D squared distance between this point a specified x, y coordinate.
A rectangle specified with double values.
QgsPointXY center() const
Returns the center point of the rectangle.
void grow(double delta)
Grows the rectangle in place by the specified amount.
A spatial index for QgsFeature objects.
QList< QgsFeatureId > intersects(const QgsRectangle &rectangle) const
Returns a list of features with a bounding box which intersects the specified rectangle.
bool addFeature(QgsFeature &feature, QgsFeatureSink::Flags flags=QgsFeatureSink::Flags()) override
Adds a feature to the index.
A class to represent a vector.
static Qgis::GeometryType geometryType(Qgis::WkbType type)
Returns the geometry type for a WKB type, e.g., both MultiPolygon and CurvePolygon would have a Polyg...
bool qgsDoubleNear(double a, double b, double epsilon=4 *std::numeric_limits< double >::epsilon())
Compare two doubles (but allow some difference)
QList< QgsFeature > QgsFeatureList
QSet< QgsFeatureId > QgsFeatureIds
qint64 QgsFeatureId
64 bit feature ids negative numbers are used for uncommitted/newly added features
void _GEOSQueryCallback(void *item, void *userdata)
QLineF segment(int index, QRectF rect, double radius)
Utility class for identifying a unique vertex within a geometry.
QList< const QgsPointCloudLayerProfileResults::PointResult * > * list