Branch data Line data Source code
1 : : /***************************************************************************
2 : : qgsgeos.h
3 : : -------------------------------------------------------------------
4 : : Date : 22 Sept 2014
5 : : Copyright : (C) 2014 by Marco Hugentobler
6 : : email : marco.hugentobler at sourcepole 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 : : #ifndef QGSGEOS_H
17 : : #define QGSGEOS_H
18 : :
19 : : #define SIP_NO_FILE
20 : :
21 : : #include "qgis_core.h"
22 : : #include "qgsgeometryengine.h"
23 : : #include "qgsgeometry.h"
24 : : #include <geos_c.h>
25 : :
26 : : #if defined(GEOS_VERSION_MAJOR) && (GEOS_VERSION_MAJOR<3)
27 : : #define GEOSGeometry struct GEOSGeom_t
28 : : #define GEOSCoordSequence struct GEOSCoordSeq_t
29 : : #endif
30 : :
31 : : class QgsLineString;
32 : : class QgsPolygon;
33 : : class QgsGeometry;
34 : : class QgsGeometryCollection;
35 : :
36 : : /**
37 : : * Contains geos related utilities and functions.
38 : : * \note not available in Python bindings.
39 : : * \since QGIS 3.0
40 : : */
41 : : namespace geos
42 : : {
43 : :
44 : : /**
45 : : * Destroys the GEOS geometry \a geom, using the static QGIS
46 : : * geos context.
47 : : */
48 : : struct GeosDeleter
49 : : {
50 : :
51 : : /**
52 : : * Destroys the GEOS geometry \a geom, using the static QGIS
53 : : * geos context.
54 : : */
55 : : void CORE_EXPORT operator()( GEOSGeometry *geom );
56 : :
57 : : /**
58 : : * Destroys the GEOS prepared geometry \a geom, using the static QGIS
59 : : * geos context.
60 : : */
61 : : void CORE_EXPORT operator()( const GEOSPreparedGeometry *geom );
62 : :
63 : : /**
64 : : * Destroys the GEOS buffer params \a params, using the static QGIS
65 : : * geos context.
66 : : */
67 : : void CORE_EXPORT operator()( GEOSBufferParams *params );
68 : :
69 : : /**
70 : : * Destroys the GEOS coordinate sequence \a sequence, using the static QGIS
71 : : * geos context.
72 : : */
73 : : void CORE_EXPORT operator()( GEOSCoordSequence *sequence );
74 : : };
75 : :
76 : : /**
77 : : * Scoped GEOS pointer.
78 : : */
79 : : using unique_ptr = std::unique_ptr< GEOSGeometry, GeosDeleter>;
80 : :
81 : : /**
82 : : * Scoped GEOS prepared geometry pointer.
83 : : */
84 : : using prepared_unique_ptr = std::unique_ptr< const GEOSPreparedGeometry, GeosDeleter>;
85 : :
86 : : /**
87 : : * Scoped GEOS buffer params pointer.
88 : : */
89 : : using buffer_params_unique_ptr = std::unique_ptr< GEOSBufferParams, GeosDeleter>;
90 : :
91 : : /**
92 : : * Scoped GEOS coordinate sequence pointer.
93 : : */
94 : : using coord_sequence_unique_ptr = std::unique_ptr< GEOSCoordSequence, GeosDeleter>;
95 : :
96 : : }
97 : :
98 : : /**
99 : : * \ingroup core
100 : : * \brief Does vector analysis using the geos library and handles import, export, exception handling*
101 : : * \note not available in Python bindings
102 : : */
103 : 30823 : class CORE_EXPORT QgsGeos: public QgsGeometryEngine
104 : : {
105 : : public:
106 : :
107 : : /**
108 : : * GEOS geometry engine constructor
109 : : * \param geometry The geometry
110 : : * \param precision The precision of the grid to which to snap the geometry vertices. If 0, no snapping is performed.
111 : : */
112 : : QgsGeos( const QgsAbstractGeometry *geometry, double precision = 0 );
113 : :
114 : : /**
115 : : * Creates a new QgsGeometry object, feeding in a geometry in GEOS format.
116 : : * This class will take ownership of the buffer.
117 : : */
118 : : static QgsGeometry geometryFromGeos( GEOSGeometry *geos );
119 : :
120 : : /**
121 : : * Creates a new QgsGeometry object, feeding in a geometry in GEOS format.
122 : : */
123 : : static QgsGeometry geometryFromGeos( const geos::unique_ptr &geos );
124 : :
125 : : /**
126 : : * Adds a new island polygon to a multipolygon feature
127 : : * \param geometry geometry to add part to
128 : : * \param newPart part to add. Ownership is NOT transferred.
129 : : * \returns OperationResult a result code: success or reason of failure
130 : : */
131 : : static QgsGeometry::OperationResult addPart( QgsGeometry &geometry, GEOSGeometry *newPart );
132 : :
133 : : void geometryChanged() override;
134 : : void prepareGeometry() override;
135 : :
136 : : QgsAbstractGeometry *intersection( const QgsAbstractGeometry *geom, QString *errorMsg = nullptr ) const override;
137 : : QgsAbstractGeometry *difference( const QgsAbstractGeometry *geom, QString *errorMsg = nullptr ) const override;
138 : :
139 : : /**
140 : : * Performs a fast, non-robust intersection between the geometry and
141 : : * a \a rectangle. The returned geometry may be invalid.
142 : : */
143 : : std::unique_ptr< QgsAbstractGeometry > clip( const QgsRectangle &rectangle, QString *errorMsg = nullptr ) const;
144 : :
145 : : /**
146 : : * Subdivides the geometry. The returned geometry will be a collection containing subdivided parts
147 : : * from the original geometry, where no part has more then the specified maximum number of nodes (\a maxNodes).
148 : : *
149 : : * This is useful for dividing a complex geometry into less complex parts, which are better able to be spatially
150 : : * indexed and faster to perform further operations such as intersects on. The returned geometry parts may
151 : : * not be valid and may contain self-intersections.
152 : : *
153 : : * The minimum allowed value for \a maxNodes is 8.
154 : : *
155 : : * Curved geometries are not supported.
156 : : *
157 : : * \since QGIS 3.0
158 : : */
159 : : std::unique_ptr< QgsAbstractGeometry > subdivide( int maxNodes, QString *errorMsg = nullptr ) const;
160 : :
161 : : QgsAbstractGeometry *combine( const QgsAbstractGeometry *geom, QString *errorMsg = nullptr ) const override;
162 : : QgsAbstractGeometry *combine( const QVector<QgsAbstractGeometry *> &geomList, QString *errorMsg ) const override;
163 : : QgsAbstractGeometry *combine( const QVector< QgsGeometry > &, QString *errorMsg = nullptr ) const override;
164 : : QgsAbstractGeometry *symDifference( const QgsAbstractGeometry *geom, QString *errorMsg = nullptr ) const override;
165 : : QgsAbstractGeometry *buffer( double distance, int segments, QString *errorMsg = nullptr ) const override;
166 : : QgsAbstractGeometry *buffer( double distance, int segments, int endCapStyle, int joinStyle, double miterLimit, QString *errorMsg = nullptr ) const override;
167 : : QgsAbstractGeometry *simplify( double tolerance, QString *errorMsg = nullptr ) const override;
168 : : QgsAbstractGeometry *interpolate( double distance, QString *errorMsg = nullptr ) const override;
169 : : QgsAbstractGeometry *envelope( QString *errorMsg = nullptr ) const override;
170 : : QgsPoint *centroid( QString *errorMsg = nullptr ) const override;
171 : : QgsPoint *pointOnSurface( QString *errorMsg = nullptr ) const override;
172 : : QgsAbstractGeometry *convexHull( QString *errorMsg = nullptr ) const override;
173 : : double distance( const QgsAbstractGeometry *geom, QString *errorMsg = nullptr ) const override;
174 : :
175 : : /**
176 : : * Returns the Hausdorff distance between this geometry and \a geom. This is basically a measure of how similar or dissimilar 2 geometries are.
177 : : *
178 : : * This algorithm is an approximation to the standard Hausdorff distance. This approximation is exact or close enough for a large
179 : : * subset of useful cases. Examples of these are:
180 : : *
181 : : * - computing distance between Linestrings that are roughly parallel to each other,
182 : : * and roughly equal in length. This occurs in matching linear networks.
183 : : * - Testing similarity of geometries.
184 : : *
185 : : * If the default approximate provided by this method is insufficient, use hausdorffDistanceDensify() instead.
186 : : *
187 : : * \see hausdorffDistanceDensify()
188 : : * \since QGIS 3.0
189 : : */
190 : : double hausdorffDistance( const QgsAbstractGeometry *geom, QString *errorMsg = nullptr ) const;
191 : :
192 : : /**
193 : : * Returns the Hausdorff distance between this geometry and \a geom. This is basically a measure of how similar or dissimilar 2 geometries are.
194 : : *
195 : : * This function accepts a \a densifyFraction argument. The function performs a segment
196 : : * densification before computing the discrete Hausdorff distance. The \a densifyFraction parameter
197 : : * sets the fraction by which to densify each segment. Each segment will be split into a
198 : : * number of equal-length subsegments, whose fraction of the total length is
199 : : * closest to the given fraction.
200 : : *
201 : : * This method can be used when the default approximation provided by hausdorffDistance()
202 : : * is not sufficient. Decreasing the \a densifyFraction parameter will make the
203 : : * distance returned approach the true Hausdorff distance for the geometries.
204 : : *
205 : : * \see hausdorffDistance()
206 : : * \since QGIS 3.0
207 : : */
208 : : double hausdorffDistanceDensify( const QgsAbstractGeometry *geom, double densifyFraction, QString *errorMsg = nullptr ) const;
209 : :
210 : : bool intersects( const QgsAbstractGeometry *geom, QString *errorMsg = nullptr ) const override;
211 : : bool touches( const QgsAbstractGeometry *geom, QString *errorMsg = nullptr ) const override;
212 : : bool crosses( const QgsAbstractGeometry *geom, QString *errorMsg = nullptr ) const override;
213 : : bool within( const QgsAbstractGeometry *geom, QString *errorMsg = nullptr ) const override;
214 : : bool overlaps( const QgsAbstractGeometry *geom, QString *errorMsg = nullptr ) const override;
215 : : bool contains( const QgsAbstractGeometry *geom, QString *errorMsg = nullptr ) const override;
216 : : bool disjoint( const QgsAbstractGeometry *geom, QString *errorMsg = nullptr ) const override;
217 : : QString relate( const QgsAbstractGeometry *geom, QString *errorMsg = nullptr ) const override;
218 : : bool relatePattern( const QgsAbstractGeometry *geom, const QString &pattern, QString *errorMsg = nullptr ) const override;
219 : : double area( QString *errorMsg = nullptr ) const override;
220 : : double length( QString *errorMsg = nullptr ) const override;
221 : : bool isValid( QString *errorMsg = nullptr, bool allowSelfTouchingHoles = false, QgsGeometry *errorLoc = nullptr ) const override;
222 : : bool isEqual( const QgsAbstractGeometry *geom, QString *errorMsg = nullptr ) const override;
223 : : bool isEmpty( QString *errorMsg = nullptr ) const override;
224 : : bool isSimple( QString *errorMsg = nullptr ) const override;
225 : :
226 : : EngineOperationResult splitGeometry( const QgsLineString &splitLine,
227 : : QVector<QgsGeometry> &newGeometries,
228 : : bool topological,
229 : : QgsPointSequence &topologyTestPoints,
230 : : QString *errorMsg = nullptr, bool skipIntersectionCheck = false ) const override;
231 : :
232 : : QgsAbstractGeometry *offsetCurve( double distance, int segments, int joinStyle, double miterLimit, QString *errorMsg = nullptr ) const override;
233 : :
234 : : /**
235 : : * Returns a single sided buffer for a geometry. The buffer is only
236 : : * applied to one side of the geometry.
237 : : * \param distance buffer distance
238 : : * \param segments for round joins, number of segments to approximate quarter-circle
239 : : * \param side side of geometry to buffer (0 = left, 1 = right)
240 : : * \param joinStyle join style for corners ( Round (1) / Miter (2) / Bevel (3) )
241 : : * \param miterLimit limit on the miter ratio used for very sharp corners
242 : : * \param errorMsg error messages emitted, if any
243 : : * \returns buffered geometry, or an NULLPTR if buffer could not be
244 : : * calculated
245 : : * \since QGIS 3.0
246 : : */
247 : : std::unique_ptr< QgsAbstractGeometry > singleSidedBuffer( double distance, int segments, int side,
248 : : int joinStyle, double miterLimit,
249 : : QString *errorMsg = nullptr ) const;
250 : :
251 : : /**
252 : : * Reshapes the geometry using a line
253 : : * \param reshapeWithLine the line used to reshape lines or polygons
254 : : * \param errorCode if specified, provides result of operation (success or reason of failure)
255 : : * \param errorMsg if specified, provides more details about failure
256 : : * \return the reshaped geometry
257 : : */
258 : : std::unique_ptr< QgsAbstractGeometry > reshapeGeometry( const QgsLineString &reshapeWithLine, EngineOperationResult *errorCode, QString *errorMsg = nullptr ) const;
259 : :
260 : : /**
261 : : * Merges any connected lines in a LineString/MultiLineString geometry and
262 : : * converts them to single line strings.
263 : : * \param errorMsg if specified, will be set to any reported GEOS errors
264 : : * \returns a LineString or MultiLineString geometry, with any connected lines
265 : : * joined. An empty geometry will be returned if the input geometry was not a
266 : : * LineString/MultiLineString geometry.
267 : : * \since QGIS 3.0
268 : : */
269 : : QgsGeometry mergeLines( QString *errorMsg = nullptr ) const;
270 : :
271 : : /**
272 : : * Returns the closest point on the geometry to the other geometry.
273 : : * \see shortestLine()
274 : : * \since QGIS 2.14
275 : : */
276 : : QgsGeometry closestPoint( const QgsGeometry &other, QString *errorMsg = nullptr ) const;
277 : :
278 : : /**
279 : : * Returns the shortest line joining this geometry to the other geometry.
280 : : * \see closestPoint()
281 : : * \since QGIS 2.14
282 : : */
283 : : QgsGeometry shortestLine( const QgsGeometry &other, QString *errorMsg = nullptr ) const;
284 : :
285 : : /**
286 : : * Returns the shortest line joining this geometry to the other geometry.
287 : : * \see closestPoint()
288 : : * \since QGIS 3.20
289 : : */
290 : : QgsGeometry shortestLine( const QgsAbstractGeometry *other, QString *errorMsg = nullptr ) const;
291 : :
292 : : /**
293 : : * Returns a distance representing the location along this linestring of the closest point
294 : : * on this linestring geometry to the specified point. Ie, the returned value indicates
295 : : * how far along this linestring you need to traverse to get to the closest location
296 : : * where this linestring comes to the specified point.
297 : : * \param point point to seek proximity to
298 : : * \param errorMsg error messages emitted, if any
299 : : * \returns distance along line, or -1 on error
300 : : * \note only valid for linestring geometries
301 : : */
302 : : double lineLocatePoint( const QgsPoint &point, QString *errorMsg = nullptr ) const;
303 : :
304 : : /**
305 : : * Creates a GeometryCollection geometry containing possible polygons formed from the constituent
306 : : * linework of a set of \a geometries. The input geometries must be fully noded (i.e. nodes exist
307 : : * at every common intersection of the geometries). The easiest way to ensure this is to first
308 : : * unary union these geometries by calling combine() on the set of input geometries and then
309 : : * pass the result to polygonize().
310 : : * An empty geometry will be returned in the case of errors.
311 : : * \since QGIS 3.0
312 : : */
313 : : static QgsGeometry polygonize( const QVector<const QgsAbstractGeometry *> &geometries, QString *errorMsg = nullptr );
314 : :
315 : : /**
316 : : * Creates a Voronoi diagram for the nodes contained within the geometry.
317 : : *
318 : : * Returns the Voronoi polygons for the nodes contained within the geometry.
319 : : * If \a extent is specified then it will be used as a clipping envelope for the diagram.
320 : : * If no extent is set then the clipping envelope will be automatically calculated.
321 : : * In either case the diagram will be clipped to the larger of the provided envelope
322 : : * OR the envelope surrounding all input nodes.
323 : : * The \a tolerance parameter specifies an optional snapping tolerance which can
324 : : * be used to improve the robustness of the diagram calculation.
325 : : * If \a edgesOnly is TRUE than line string boundary geometries will be returned
326 : : * instead of polygons.
327 : : * An empty geometry will be returned if the diagram could not be calculated.
328 : : * \since QGIS 3.0
329 : : */
330 : : QgsGeometry voronoiDiagram( const QgsAbstractGeometry *extent = nullptr, double tolerance = 0.0, bool edgesOnly = false, QString *errorMsg = nullptr ) const;
331 : :
332 : : /**
333 : : * Returns the Delaunay triangulation for the vertices of the geometry.
334 : : * The \a tolerance parameter specifies an optional snapping tolerance which can
335 : : * be used to improve the robustness of the triangulation.
336 : : * If \a edgesOnly is TRUE than line string boundary geometries will be returned
337 : : * instead of polygons.
338 : : * An empty geometry will be returned if the diagram could not be calculated.
339 : : * \since QGIS 3.0
340 : : */
341 : : QgsGeometry delaunayTriangulation( double tolerance = 0.0, bool edgesOnly = false, QString *errorMsg = nullptr ) const;
342 : :
343 : : /**
344 : : * Create a geometry from a GEOSGeometry
345 : : * \param geos GEOSGeometry. Ownership is NOT transferred.
346 : : */
347 : : static std::unique_ptr< QgsAbstractGeometry > fromGeos( const GEOSGeometry *geos );
348 : : static std::unique_ptr< QgsPolygon > fromGeosPolygon( const GEOSGeometry *geos );
349 : :
350 : :
351 : : /**
352 : : * Returns a geos geometry - caller takes ownership of the object (should be deleted with GEOSGeom_destroy_r)
353 : : * \param geometry geometry to convert to GEOS representation
354 : : * \param precision The precision of the grid to which to snap the geometry vertices. If 0, no snapping is performed.
355 : : */
356 : : static geos::unique_ptr asGeos( const QgsGeometry &geometry, double precision = 0 );
357 : :
358 : : /**
359 : : * Returns a geos geometry - caller takes ownership of the object (should be deleted with GEOSGeom_destroy_r)
360 : : * \param geometry geometry to convert to GEOS representation
361 : : * \param precision The precision of the grid to which to snap the geometry vertices. If 0, no snapping is performed.
362 : : */
363 : : static geos::unique_ptr asGeos( const QgsAbstractGeometry *geometry, double precision = 0 );
364 : : static QgsPoint coordSeqPoint( const GEOSCoordSequence *cs, int i, bool hasZ, bool hasM );
365 : :
366 : : static GEOSContextHandle_t getGEOSHandler();
367 : :
368 : :
369 : : private:
370 : : mutable geos::unique_ptr mGeos;
371 : : geos::prepared_unique_ptr mGeosPrepared;
372 : : double mPrecision = 0.0;
373 : :
374 : : enum Overlay
375 : : {
376 : : OverlayIntersection,
377 : : OverlayDifference,
378 : : OverlayUnion,
379 : : OverlaySymDifference
380 : : };
381 : :
382 : : enum Relation
383 : : {
384 : : RelationIntersects,
385 : : RelationTouches,
386 : : RelationCrosses,
387 : : RelationWithin,
388 : : RelationOverlaps,
389 : : RelationContains,
390 : : RelationDisjoint
391 : : };
392 : :
393 : : //geos util functions
394 : : void cacheGeos() const;
395 : : std::unique_ptr< QgsAbstractGeometry > overlay( const QgsAbstractGeometry *geom, Overlay op, QString *errorMsg = nullptr ) const;
396 : : bool relation( const QgsAbstractGeometry *geom, Relation r, QString *errorMsg = nullptr ) const;
397 : : static GEOSCoordSequence *createCoordinateSequence( const QgsCurve *curve, double precision, bool forceClose = false );
398 : : static std::unique_ptr< QgsLineString > sequenceToLinestring( const GEOSGeometry *geos, bool hasZ, bool hasM );
399 : : static int numberOfGeometries( GEOSGeometry *g );
400 : : static geos::unique_ptr nodeGeometries( const GEOSGeometry *splitLine, const GEOSGeometry *geom );
401 : : int mergeGeometriesMultiTypeSplit( QVector<GEOSGeometry *> &splitResult ) const;
402 : :
403 : : /**
404 : : * Ownership of geoms is transferred
405 : : */
406 : : static geos::unique_ptr createGeosCollection( int typeId, const QVector<GEOSGeometry *> &geoms );
407 : :
408 : : static geos::unique_ptr createGeosPointXY( double x, double y, bool hasZ, double z, bool hasM, double m, int coordDims, double precision );
409 : : static geos::unique_ptr createGeosPoint( const QgsAbstractGeometry *point, int coordDims, double precision );
410 : : static geos::unique_ptr createGeosLinestring( const QgsAbstractGeometry *curve, double precision );
411 : : static geos::unique_ptr createGeosPolygon( const QgsAbstractGeometry *poly, double precision );
412 : :
413 : : //utils for geometry split
414 : : bool topologicalTestPointsSplit( const GEOSGeometry *splitLine, QgsPointSequence &testPoints, QString *errorMsg = nullptr ) const;
415 : : geos::unique_ptr linePointDifference( GEOSGeometry *GEOSsplitPoint ) const;
416 : : EngineOperationResult splitLinearGeometry( GEOSGeometry *splitLine, QVector<QgsGeometry > &newGeometries, bool skipIntersectionCheck ) const;
417 : : EngineOperationResult splitPolygonGeometry( GEOSGeometry *splitLine, QVector<QgsGeometry > &newGeometries, bool skipIntersectionCheck ) const;
418 : :
419 : : //utils for reshape
420 : : static geos::unique_ptr reshapeLine( const GEOSGeometry *line, const GEOSGeometry *reshapeLineGeos, double precision );
421 : : static geos::unique_ptr reshapePolygon( const GEOSGeometry *polygon, const GEOSGeometry *reshapeLineGeos, double precision );
422 : : static int lineContainedInLine( const GEOSGeometry *line1, const GEOSGeometry *line2 );
423 : : static int pointContainedInLine( const GEOSGeometry *point, const GEOSGeometry *line );
424 : : static int geomDigits( const GEOSGeometry *geom );
425 : : void subdivideRecursive( const GEOSGeometry *currentPart, int maxNodes, int depth, QgsGeometryCollection *parts, const QgsRectangle &clipRect ) const;
426 : : };
427 : :
428 : : /// @cond PRIVATE
429 : :
430 : :
431 : 9 : class GEOSException : public std::runtime_error
432 : : {
433 : : public:
434 : 9 : explicit GEOSException( const QString &message )
435 : 9 : : std::runtime_error( message.toUtf8().constData() )
436 : 18 : {
437 : 9 : }
438 : : };
439 : :
440 : : /// @endcond
441 : :
442 : : #endif // QGSGEOS_H
|