QGIS API Documentation  2.14.0-Essen
qgsgeos.cpp
Go to the documentation of this file.
1 /***************************************************************************
2  qgsgeos.cpp
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 #include "qgsgeos.h"
17 #include "qgsabstractgeometryv2.h"
19 #include "qgsgeometryfactory.h"
20 #include "qgslinestringv2.h"
21 #include "qgsmessagelog.h"
22 #include "qgsmulticurvev2.h"
23 #include "qgsmultilinestringv2.h"
24 #include "qgsmultipointv2.h"
25 #include "qgsmultipolygonv2.h"
26 #include "qgslogger.h"
27 #include "qgspolygonv2.h"
28 #include <limits>
29 #include <cstdio>
30 #include <QtCore/qmath.h>
31 
32 #define DEFAULT_QUADRANT_SEGMENTS 8
33 
34 #define CATCH_GEOS(r) \
35  catch (GEOSException &e) \
36  { \
37  QgsMessageLog::logMessage( QObject::tr( "Exception: %1" ).arg( e.what() ), QObject::tr("GEOS") ); \
38  return r; \
39  }
40 
41 #define CATCH_GEOS_WITH_ERRMSG(r) \
42  catch (GEOSException &e) \
43  { \
44  QgsMessageLog::logMessage( QObject::tr( "Exception: %1" ).arg( e.what() ), QObject::tr("GEOS") ); \
45  if ( errorMsg ) \
46  { \
47  *errorMsg = e.what(); \
48  } \
49  return r; \
50  }
51 
53 
54 static void throwGEOSException( const char *fmt, ... )
55 {
56  va_list ap;
57  char buffer[1024];
58 
59  va_start( ap, fmt );
60  vsnprintf( buffer, sizeof buffer, fmt, ap );
61  va_end( ap );
62 
63  QgsDebugMsg( QString( "GEOS exception: %1" ).arg( buffer ) );
64 
65  throw GEOSException( QString::fromUtf8( buffer ) );
66 }
67 
68 static void printGEOSNotice( const char *fmt, ... )
69 {
70 #if defined(QGISDEBUG)
71  va_list ap;
72  char buffer[1024];
73 
74  va_start( ap, fmt );
75  vsnprintf( buffer, sizeof buffer, fmt, ap );
76  va_end( ap );
77 
78  QgsDebugMsg( QString( "GEOS notice: %1" ).arg( QString::fromUtf8( buffer ) ) );
79 #else
80  Q_UNUSED( fmt );
81 #endif
82 }
83 
84 class GEOSInit
85 {
86  public:
87  GEOSContextHandle_t ctxt;
88 
89  GEOSInit()
90  {
91  ctxt = initGEOS_r( printGEOSNotice, throwGEOSException );
92  }
93 
94  ~GEOSInit()
95  {
96  finishGEOS_r( ctxt );
97  }
98 
99  private:
100 
101  GEOSInit( const GEOSInit& rh );
102  GEOSInit& operator=( const GEOSInit& rh );
103 };
104 
105 static GEOSInit geosinit;
106 
108 
109 
115 {
116  public:
117  explicit GEOSGeomScopedPtr( GEOSGeometry* geom = nullptr ) : mGeom( geom ) {}
118  ~GEOSGeomScopedPtr() { GEOSGeom_destroy_r( geosinit.ctxt, mGeom ); }
119  GEOSGeometry* get() const { return mGeom; }
120  operator bool() const { return nullptr != mGeom; }
121  void reset( GEOSGeometry* geom )
122  {
123  GEOSGeom_destroy_r( geosinit.ctxt, mGeom );
124  mGeom = geom;
125  }
126 
127  private:
128  GEOSGeometry* mGeom;
129 
130  private:
132  GEOSGeomScopedPtr& operator=( const GEOSGeomScopedPtr& rh );
133 };
134 
135 QgsGeos::QgsGeos( const QgsAbstractGeometryV2* geometry, double precision )
136  : QgsGeometryEngine( geometry ), mGeos( nullptr ), mGeosPrepared( nullptr ), mPrecision( precision )
137 {
138  cacheGeos();
139 }
140 
142 {
143  GEOSGeom_destroy_r( geosinit.ctxt, mGeos );
144  mGeos = nullptr;
145  GEOSPreparedGeom_destroy_r( geosinit.ctxt, mGeosPrepared );
146  mGeosPrepared = nullptr;
147 }
148 
150 {
151  GEOSGeom_destroy_r( geosinit.ctxt, mGeos );
152  mGeos = nullptr;
153  GEOSPreparedGeom_destroy_r( geosinit.ctxt, mGeosPrepared );
154  mGeosPrepared = nullptr;
155  cacheGeos();
156 }
157 
159 {
160  GEOSPreparedGeom_destroy_r( geosinit.ctxt, mGeosPrepared );
161  mGeosPrepared = nullptr;
162  if ( mGeos )
163  {
164  mGeosPrepared = GEOSPrepare_r( geosinit.ctxt, mGeos );
165  }
166 }
167 
168 void QgsGeos::cacheGeos() const
169 {
170  if ( !mGeometry || mGeos )
171  {
172  return;
173  }
174 
175  mGeos = asGeos( mGeometry, mPrecision );
176 }
177 
179 {
180  return overlay( geom, INTERSECTION, errorMsg );
181 }
182 
184 {
185  return overlay( geom, DIFFERENCE, errorMsg );
186 }
187 
189 {
190  return overlay( geom, UNION, errorMsg );
191 }
192 
194 {
195 
196  QVector< GEOSGeometry* > geosGeometries;
197  geosGeometries.resize( geomList.size() );
198  for ( int i = 0; i < geomList.size(); ++i )
199  {
200  geosGeometries[i] = asGeos( geomList.at( i ), mPrecision );
201  }
202 
203  GEOSGeometry* geomUnion = nullptr;
204  try
205  {
206  GEOSGeometry* geomCollection = createGeosCollection( GEOS_GEOMETRYCOLLECTION, geosGeometries );
207  geomUnion = GEOSUnaryUnion_r( geosinit.ctxt, geomCollection );
208  GEOSGeom_destroy_r( geosinit.ctxt, geomCollection );
209  }
210  CATCH_GEOS_WITH_ERRMSG( nullptr )
211 
212  QgsAbstractGeometryV2* result = fromGeos( geomUnion );
213  GEOSGeom_destroy_r( geosinit.ctxt, geomUnion );
214  return result;
215 }
216 
218 {
219  return overlay( geom, SYMDIFFERENCE, errorMsg );
220 }
221 
222 double QgsGeos::distance( const QgsAbstractGeometryV2& geom, QString* errorMsg ) const
223 {
224  double distance = -1.0;
225  if ( !mGeos )
226  {
227  return distance;
228  }
229 
230  GEOSGeometry* otherGeosGeom = asGeos( &geom, mPrecision );
231  if ( !otherGeosGeom )
232  {
233  return distance;
234  }
235 
236  try
237  {
238  GEOSDistance_r( geosinit.ctxt, mGeos, otherGeosGeom, &distance );
239  }
240  CATCH_GEOS_WITH_ERRMSG( -1.0 )
241 
242  GEOSGeom_destroy_r( geosinit.ctxt, otherGeosGeom );
243 
244  return distance;
245 }
246 
247 bool QgsGeos::intersects( const QgsAbstractGeometryV2& geom, QString* errorMsg ) const
248 {
249  return relation( geom, INTERSECTS, errorMsg );
250 }
251 
252 bool QgsGeos::touches( const QgsAbstractGeometryV2& geom, QString* errorMsg ) const
253 {
254  return relation( geom, TOUCHES, errorMsg );
255 }
256 
257 bool QgsGeos::crosses( const QgsAbstractGeometryV2& geom, QString* errorMsg ) const
258 {
259  return relation( geom, CROSSES, errorMsg );
260 }
261 
262 bool QgsGeos::within( const QgsAbstractGeometryV2& geom, QString* errorMsg ) const
263 {
264  return relation( geom, WITHIN, errorMsg );
265 }
266 
267 bool QgsGeos::overlaps( const QgsAbstractGeometryV2& geom, QString* errorMsg ) const
268 {
269  return relation( geom, OVERLAPS, errorMsg );
270 }
271 
272 bool QgsGeos::contains( const QgsAbstractGeometryV2& geom, QString* errorMsg ) const
273 {
274  return relation( geom, CONTAINS, errorMsg );
275 }
276 
277 bool QgsGeos::disjoint( const QgsAbstractGeometryV2& geom, QString* errorMsg ) const
278 {
279  return relation( geom, DISJOINT, errorMsg );
280 }
281 
282 QString QgsGeos::relate( const QgsAbstractGeometryV2& geom, QString* errorMsg ) const
283 {
284  if ( !mGeos )
285  {
286  return QString();
287  }
288 
289  GEOSGeomScopedPtr geosGeom( asGeos( &geom, mPrecision ) );
290  if ( !geosGeom )
291  {
292  return QString();
293  }
294 
295  QString result;
296  try
297  {
298  char* r = GEOSRelate_r( geosinit.ctxt, mGeos, geosGeom.get() );
299  if ( r )
300  {
301  result = QString( r );
302  GEOSFree_r( geosinit.ctxt, r );
303  }
304  }
305  catch ( GEOSException &e )
306  {
307  if ( errorMsg )
308  {
309  *errorMsg = e.what();
310  }
311  }
312 
313  return result;
314 }
315 
316 bool QgsGeos::relatePattern( const QgsAbstractGeometryV2& geom, const QString& pattern, QString* errorMsg ) const
317 {
318  if ( !mGeos )
319  {
320  return false;
321  }
322 
323  GEOSGeomScopedPtr geosGeom( asGeos( &geom, mPrecision ) );
324  if ( !geosGeom )
325  {
326  return false;
327  }
328 
329  bool result = false;
330  try
331  {
332  result = ( GEOSRelatePattern_r( geosinit.ctxt, mGeos, geosGeom.get(), pattern.toLocal8Bit().constData() ) == 1 );
333  }
334  catch ( GEOSException &e )
335  {
336  if ( errorMsg )
337  {
338  *errorMsg = e.what();
339  }
340  }
341 
342  return result;
343 }
344 
345 double QgsGeos::area( QString* errorMsg ) const
346 {
347  double area = -1.0;
348  if ( !mGeos )
349  {
350  return area;
351  }
352 
353  try
354  {
355  if ( GEOSArea_r( geosinit.ctxt, mGeos, &area ) != 1 )
356  return -1.0;
357  }
358  CATCH_GEOS_WITH_ERRMSG( -1.0 );
359  return area;
360 }
361 
362 double QgsGeos::length( QString* errorMsg ) const
363 {
364  double length = -1.0;
365  if ( !mGeos )
366  {
367  return length;
368  }
369  try
370  {
371  if ( GEOSLength_r( geosinit.ctxt, mGeos, &length ) != 1 )
372  return -1.0;
373  }
374  CATCH_GEOS_WITH_ERRMSG( -1.0 )
375  return length;
376 }
377 
379  QList<QgsAbstractGeometryV2*>&newGeometries,
380  bool topological,
381  QgsPointSequenceV2 &topologyTestPoints,
382  QString* errorMsg ) const
383 {
384 
385  int returnCode = 0;
386  if ( !mGeometry || !mGeos )
387  {
388  return 1;
389  }
390 
391  //return if this type is point/multipoint
392  if ( mGeometry->dimension() == 0 )
393  {
394  return 1; //cannot split points
395  }
396 
397  if ( !GEOSisValid_r( geosinit.ctxt, mGeos ) )
398  return 7;
399 
400  //make sure splitLine is valid
401  if (( mGeometry->dimension() == 1 && splitLine.numPoints() < 1 ) ||
402  ( mGeometry->dimension() == 2 && splitLine.numPoints() < 2 ) )
403  return 1;
404 
405  newGeometries.clear();
406  GEOSGeometry* splitLineGeos = nullptr;
407 
408  try
409  {
410  if ( splitLine.numPoints() > 1 )
411  {
412  splitLineGeos = createGeosLinestring( &splitLine, mPrecision );
413  }
414  else if ( splitLine.numPoints() == 1 )
415  {
416  QgsPointV2 pt = splitLine.pointN( 0 );
417  splitLineGeos = createGeosPoint( &pt, 2, mPrecision );
418  }
419  else
420  {
421  return 1;
422  }
423 
424  if ( !GEOSisValid_r( geosinit.ctxt, splitLineGeos ) || !GEOSisSimple_r( geosinit.ctxt, splitLineGeos ) )
425  {
426  GEOSGeom_destroy_r( geosinit.ctxt, splitLineGeos );
427  return 1;
428  }
429 
430  if ( topological )
431  {
432  //find out candidate points for topological corrections
433  if ( topologicalTestPointsSplit( splitLineGeos, topologyTestPoints ) != 0 )
434  return 1;
435  }
436 
437  //call split function depending on geometry type
438  if ( mGeometry->dimension() == 1 )
439  {
440  returnCode = splitLinearGeometry( splitLineGeos, newGeometries );
441  GEOSGeom_destroy_r( geosinit.ctxt, splitLineGeos );
442  }
443  else if ( mGeometry->dimension() == 2 )
444  {
445  returnCode = splitPolygonGeometry( splitLineGeos, newGeometries );
446  GEOSGeom_destroy_r( geosinit.ctxt, splitLineGeos );
447  }
448  else
449  {
450  return 1;
451  }
452  }
454 
455  return returnCode;
456 }
457 
458 
459 
460 int QgsGeos::topologicalTestPointsSplit( const GEOSGeometry* splitLine, QgsPointSequenceV2 &testPoints, QString* errorMsg ) const
461 {
462  //Find out the intersection points between splitLineGeos and this geometry.
463  //These points need to be tested for topological correctness by the calling function
464  //if topological editing is enabled
465 
466  if ( !mGeos )
467  {
468  return 1;
469  }
470 
471  try
472  {
473  testPoints.clear();
474  GEOSGeometry* intersectionGeom = GEOSIntersection_r( geosinit.ctxt, mGeos, splitLine );
475  if ( !intersectionGeom )
476  return 1;
477 
478  bool simple = false;
479  int nIntersectGeoms = 1;
480  if ( GEOSGeomTypeId_r( geosinit.ctxt, intersectionGeom ) == GEOS_LINESTRING
481  || GEOSGeomTypeId_r( geosinit.ctxt, intersectionGeom ) == GEOS_POINT )
482  simple = true;
483 
484  if ( !simple )
485  nIntersectGeoms = GEOSGetNumGeometries_r( geosinit.ctxt, intersectionGeom );
486 
487  for ( int i = 0; i < nIntersectGeoms; ++i )
488  {
489  const GEOSGeometry* currentIntersectGeom;
490  if ( simple )
491  currentIntersectGeom = intersectionGeom;
492  else
493  currentIntersectGeom = GEOSGetGeometryN_r( geosinit.ctxt, intersectionGeom, i );
494 
495  const GEOSCoordSequence* lineSequence = GEOSGeom_getCoordSeq_r( geosinit.ctxt, currentIntersectGeom );
496  unsigned int sequenceSize = 0;
497  double x, y;
498  if ( GEOSCoordSeq_getSize_r( geosinit.ctxt, lineSequence, &sequenceSize ) != 0 )
499  {
500  for ( unsigned int i = 0; i < sequenceSize; ++i )
501  {
502  if ( GEOSCoordSeq_getX_r( geosinit.ctxt, lineSequence, i, &x ) != 0 )
503  {
504  if ( GEOSCoordSeq_getY_r( geosinit.ctxt, lineSequence, i, &y ) != 0 )
505  {
506  testPoints.push_back( QgsPointV2( x, y ) );
507  }
508  }
509  }
510  }
511  }
512  GEOSGeom_destroy_r( geosinit.ctxt, intersectionGeom );
513  }
515 
516  return 0;
517 }
518 
519 GEOSGeometry* QgsGeos::linePointDifference( GEOSGeometry* GEOSsplitPoint ) const
520 {
521  int type = GEOSGeomTypeId_r( geosinit.ctxt, mGeos );
522 
523  QgsMultiCurveV2* multiCurve = nullptr;
524  if ( type == GEOS_MULTILINESTRING )
525  {
526  multiCurve = dynamic_cast<QgsMultiCurveV2*>( mGeometry->clone() );
527  }
528  else if ( type == GEOS_LINESTRING )
529  {
530  multiCurve = new QgsMultiCurveV2();
531  multiCurve->addGeometry( mGeometry->clone() );
532  }
533  else
534  {
535  return nullptr;
536  }
537 
538  if ( !multiCurve )
539  {
540  return nullptr;
541  }
542 
543 
544  QgsAbstractGeometryV2* splitGeom = fromGeos( GEOSsplitPoint );
545  QgsPointV2* splitPoint = dynamic_cast<QgsPointV2*>( splitGeom );
546  if ( !splitPoint )
547  {
548  delete splitGeom;
549  return nullptr;
550  }
551 
552  QgsMultiCurveV2 lines;
553 
554  //For each part
555  for ( int i = 0; i < multiCurve->numGeometries(); ++i )
556  {
557  const QgsLineStringV2* line = dynamic_cast<const QgsLineStringV2*>( multiCurve->geometryN( i ) );
558  if ( line )
559  {
560  //For each segment
561  QgsLineStringV2 newLine;
562  newLine.addVertex( line->pointN( 0 ) );
563  int nVertices = line->numPoints();
564  for ( int j = 1; j < ( nVertices - 1 ); ++j )
565  {
566  QgsPointV2 currentPoint = line->pointN( j );
567  newLine.addVertex( currentPoint );
568  if ( currentPoint == *splitPoint )
569  {
570  lines.addGeometry( newLine.clone() );
571  newLine = QgsLineStringV2();
572  newLine.addVertex( currentPoint );
573  }
574  }
575  newLine.addVertex( line->pointN( nVertices - 1 ) );
576  lines.addGeometry( newLine.clone() );
577  }
578  }
579 
580  delete splitGeom;
581  delete multiCurve;
582  return asGeos( &lines, mPrecision );
583 }
584 
585 int QgsGeos::splitLinearGeometry( GEOSGeometry* splitLine, QList<QgsAbstractGeometryV2*>& newGeometries ) const
586 {
587  if ( !splitLine )
588  return 2;
589 
590  if ( !mGeos )
591  return 5;
592 
593  //first test if linestring intersects geometry. If not, return straight away
594  if ( !GEOSIntersects_r( geosinit.ctxt, splitLine, mGeos ) )
595  return 1;
596 
597  //check that split line has no linear intersection
598  int linearIntersect = GEOSRelatePattern_r( geosinit.ctxt, mGeos, splitLine, "1********" );
599  if ( linearIntersect > 0 )
600  return 3;
601 
602  int splitGeomType = GEOSGeomTypeId_r( geosinit.ctxt, splitLine );
603 
604  GEOSGeometry* splitGeom;
605  if ( splitGeomType == GEOS_POINT )
606  {
607  splitGeom = linePointDifference( splitLine );
608  }
609  else
610  {
611  splitGeom = GEOSDifference_r( geosinit.ctxt, mGeos, splitLine );
612  }
613  QVector<GEOSGeometry*> lineGeoms;
614 
615  int splitType = GEOSGeomTypeId_r( geosinit.ctxt, splitGeom );
616  if ( splitType == GEOS_MULTILINESTRING )
617  {
618  int nGeoms = GEOSGetNumGeometries_r( geosinit.ctxt, splitGeom );
619  lineGeoms.reserve( nGeoms );
620  for ( int i = 0; i < nGeoms; ++i )
621  lineGeoms << GEOSGeom_clone_r( geosinit.ctxt, GEOSGetGeometryN_r( geosinit.ctxt, splitGeom, i ) );
622 
623  }
624  else
625  {
626  lineGeoms << GEOSGeom_clone_r( geosinit.ctxt, splitGeom );
627  }
628 
629  mergeGeometriesMultiTypeSplit( lineGeoms );
630 
631  for ( int i = 0; i < lineGeoms.size(); ++i )
632  {
633  newGeometries << fromGeos( lineGeoms[i] );
634  GEOSGeom_destroy_r( geosinit.ctxt, lineGeoms[i] );
635  }
636 
637  GEOSGeom_destroy_r( geosinit.ctxt, splitGeom );
638  return 0;
639 }
640 
641 int QgsGeos::splitPolygonGeometry( GEOSGeometry* splitLine, QList<QgsAbstractGeometryV2*>& newGeometries ) const
642 {
643  if ( !splitLine )
644  return 2;
645 
646  if ( !mGeos )
647  return 5;
648 
649  //first test if linestring intersects geometry. If not, return straight away
650  if ( !GEOSIntersects_r( geosinit.ctxt, splitLine, mGeos ) )
651  return 1;
652 
653  //first union all the polygon rings together (to get them noded, see JTS developer guide)
654  GEOSGeometry *nodedGeometry = nodeGeometries( splitLine, mGeos );
655  if ( !nodedGeometry )
656  return 2; //an error occurred during noding
657 
658  GEOSGeometry *polygons = GEOSPolygonize_r( geosinit.ctxt, &nodedGeometry, 1 );
659  if ( !polygons || numberOfGeometries( polygons ) == 0 )
660  {
661  if ( polygons )
662  GEOSGeom_destroy_r( geosinit.ctxt, polygons );
663 
664  GEOSGeom_destroy_r( geosinit.ctxt, nodedGeometry );
665 
666  return 4;
667  }
668 
669  GEOSGeom_destroy_r( geosinit.ctxt, nodedGeometry );
670 
671  //test every polygon if contained in original geometry
672  //include in result if yes
673  QVector<GEOSGeometry*> testedGeometries;
674  GEOSGeometry *intersectGeometry = nullptr;
675 
676  //ratio intersect geometry / geometry. This should be close to 1
677  //if the polygon belongs to the input geometry
678 
679  for ( int i = 0; i < numberOfGeometries( polygons ); i++ )
680  {
681  const GEOSGeometry *polygon = GEOSGetGeometryN_r( geosinit.ctxt, polygons, i );
682  intersectGeometry = GEOSIntersection_r( geosinit.ctxt, mGeos, polygon );
683  if ( !intersectGeometry )
684  {
685  QgsDebugMsg( "intersectGeometry is nullptr" );
686  continue;
687  }
688 
689  double intersectionArea;
690  GEOSArea_r( geosinit.ctxt, intersectGeometry, &intersectionArea );
691 
692  double polygonArea;
693  GEOSArea_r( geosinit.ctxt, polygon, &polygonArea );
694 
695  const double areaRatio = intersectionArea / polygonArea;
696  if ( areaRatio > 0.99 && areaRatio < 1.01 )
697  testedGeometries << GEOSGeom_clone_r( geosinit.ctxt, polygon );
698 
699  GEOSGeom_destroy_r( geosinit.ctxt, intersectGeometry );
700  }
701  GEOSGeom_destroy_r( geosinit.ctxt, polygons );
702 
703  bool splitDone = true;
704  int nGeometriesThis = numberOfGeometries( mGeos ); //original number of geometries
705  if ( testedGeometries.size() == nGeometriesThis )
706  {
707  splitDone = false;
708  }
709 
710  mergeGeometriesMultiTypeSplit( testedGeometries );
711 
712  //no split done, preserve original geometry
713  if ( !splitDone )
714  {
715  for ( int i = 0; i < testedGeometries.size(); ++i )
716  {
717  GEOSGeom_destroy_r( geosinit.ctxt, testedGeometries[i] );
718  }
719  return 1;
720  }
721 
722  int i;
723  for ( i = 0; i < testedGeometries.size() && GEOSisValid_r( geosinit.ctxt, testedGeometries[i] ); ++i )
724  ;
725 
726  if ( i < testedGeometries.size() )
727  {
728  for ( i = 0; i < testedGeometries.size(); ++i )
729  GEOSGeom_destroy_r( geosinit.ctxt, testedGeometries[i] );
730 
731  return 3;
732  }
733 
734  for ( i = 0; i < testedGeometries.size(); ++i )
735  newGeometries << fromGeos( testedGeometries[i] );
736 
737  return 0;
738 }
739 
740 GEOSGeometry* QgsGeos::nodeGeometries( const GEOSGeometry *splitLine, const GEOSGeometry *geom )
741 {
742  if ( !splitLine || !geom )
743  return nullptr;
744 
745  GEOSGeometry *geometryBoundary = nullptr;
746  if ( GEOSGeomTypeId_r( geosinit.ctxt, geom ) == GEOS_POLYGON || GEOSGeomTypeId_r( geosinit.ctxt, geom ) == GEOS_MULTIPOLYGON )
747  geometryBoundary = GEOSBoundary_r( geosinit.ctxt, geom );
748  else
749  geometryBoundary = GEOSGeom_clone_r( geosinit.ctxt, geom );
750 
751  GEOSGeometry *splitLineClone = GEOSGeom_clone_r( geosinit.ctxt, splitLine );
752  GEOSGeometry *unionGeometry = GEOSUnion_r( geosinit.ctxt, splitLineClone, geometryBoundary );
753  GEOSGeom_destroy_r( geosinit.ctxt, splitLineClone );
754 
755  GEOSGeom_destroy_r( geosinit.ctxt, geometryBoundary );
756  return unionGeometry;
757 }
758 
759 int QgsGeos::mergeGeometriesMultiTypeSplit( QVector<GEOSGeometry*>& splitResult ) const
760 {
761  if ( !mGeos )
762  return 1;
763 
764  //convert mGeos to geometry collection
765  int type = GEOSGeomTypeId_r( geosinit.ctxt, mGeos );
766  if ( type != GEOS_GEOMETRYCOLLECTION &&
767  type != GEOS_MULTILINESTRING &&
768  type != GEOS_MULTIPOLYGON &&
769  type != GEOS_MULTIPOINT )
770  return 0;
771 
772  QVector<GEOSGeometry*> copyList = splitResult;
773  splitResult.clear();
774 
775  //collect all the geometries that belong to the initial multifeature
776  QVector<GEOSGeometry*> unionGeom;
777 
778  for ( int i = 0; i < copyList.size(); ++i )
779  {
780  //is this geometry a part of the original multitype?
781  bool isPart = false;
782  for ( int j = 0; j < GEOSGetNumGeometries_r( geosinit.ctxt, mGeos ); j++ )
783  {
784  if ( GEOSEquals_r( geosinit.ctxt, copyList[i], GEOSGetGeometryN_r( geosinit.ctxt, mGeos, j ) ) )
785  {
786  isPart = true;
787  break;
788  }
789  }
790 
791  if ( isPart )
792  {
793  unionGeom << copyList[i];
794  }
795  else
796  {
797  QVector<GEOSGeometry*> geomVector;
798  geomVector << copyList[i];
799 
800  if ( type == GEOS_MULTILINESTRING )
801  splitResult << createGeosCollection( GEOS_MULTILINESTRING, geomVector );
802  else if ( type == GEOS_MULTIPOLYGON )
803  splitResult << createGeosCollection( GEOS_MULTIPOLYGON, geomVector );
804  else
805  GEOSGeom_destroy_r( geosinit.ctxt, copyList[i] );
806  }
807  }
808 
809  //make multifeature out of unionGeom
810  if ( !unionGeom.isEmpty() )
811  {
812  if ( type == GEOS_MULTILINESTRING )
813  splitResult << createGeosCollection( GEOS_MULTILINESTRING, unionGeom );
814  else if ( type == GEOS_MULTIPOLYGON )
815  splitResult << createGeosCollection( GEOS_MULTIPOLYGON, unionGeom );
816  }
817  else
818  {
819  unionGeom.clear();
820  }
821 
822  return 0;
823 }
824 
825 GEOSGeometry* QgsGeos::createGeosCollection( int typeId, const QVector<GEOSGeometry*>& geoms )
826 {
827  int nNullGeoms = geoms.count( nullptr );
828  int nNotNullGeoms = geoms.size() - nNullGeoms;
829 
830  GEOSGeometry **geomarr = new GEOSGeometry*[ nNotNullGeoms ];
831  if ( !geomarr )
832  {
833  return nullptr;
834  }
835 
836  int i = 0;
838  for ( ; geomIt != geoms.constEnd(); ++geomIt )
839  {
840  if ( *geomIt )
841  {
842  geomarr[i] = *geomIt;
843  ++i;
844  }
845  }
846  GEOSGeometry *geom = nullptr;
847 
848  try
849  {
850  geom = GEOSGeom_createCollection_r( geosinit.ctxt, typeId, geomarr, nNotNullGeoms );
851  }
852  catch ( GEOSException &e )
853  {
854  QgsMessageLog::logMessage( QObject::tr( "Exception: %1" ).arg( e.what() ), QObject::tr( "GEOS" ) );
855  }
856 
857  delete [] geomarr;
858 
859  return geom;
860 }
861 
862 QgsAbstractGeometryV2* QgsGeos::fromGeos( const GEOSGeometry* geos )
863 {
864  if ( !geos )
865  {
866  return nullptr;
867  }
868 
869  int nCoordDims = GEOSGeom_getCoordinateDimension_r( geosinit.ctxt, geos );
870  int nDims = GEOSGeom_getDimensions_r( geosinit.ctxt, geos );
871  bool hasZ = ( nCoordDims == 3 );
872  bool hasM = (( nDims - nCoordDims ) == 1 );
873 
874  switch ( GEOSGeomTypeId_r( geosinit.ctxt, geos ) )
875  {
876  case GEOS_POINT: // a point
877  {
878  const GEOSCoordSequence* cs = GEOSGeom_getCoordSeq_r( geosinit.ctxt, geos );
879  return ( coordSeqPoint( cs, 0, hasZ, hasM ).clone() );
880  }
881  case GEOS_LINESTRING:
882  {
883  return sequenceToLinestring( geos, hasZ, hasM );
884  }
885  case GEOS_POLYGON:
886  {
887  return fromGeosPolygon( geos );
888  }
889  case GEOS_MULTIPOINT:
890  {
891  QgsMultiPointV2* multiPoint = new QgsMultiPointV2();
892  int nParts = GEOSGetNumGeometries_r( geosinit.ctxt, geos );
893  for ( int i = 0; i < nParts; ++i )
894  {
895  const GEOSCoordSequence* cs = GEOSGeom_getCoordSeq_r( geosinit.ctxt, GEOSGetGeometryN_r( geosinit.ctxt, geos, i ) );
896  if ( cs )
897  {
898  multiPoint->addGeometry( coordSeqPoint( cs, 0, hasZ, hasM ).clone() );
899  }
900  }
901  return multiPoint;
902  }
903  case GEOS_MULTILINESTRING:
904  {
905  QgsMultiLineStringV2* multiLineString = new QgsMultiLineStringV2();
906  int nParts = GEOSGetNumGeometries_r( geosinit.ctxt, geos );
907  for ( int i = 0; i < nParts; ++i )
908  {
909  QgsLineStringV2* line = sequenceToLinestring( GEOSGetGeometryN_r( geosinit.ctxt, geos, i ), hasZ, hasM );
910  if ( line )
911  {
912  multiLineString->addGeometry( line );
913  }
914  }
915  return multiLineString;
916  }
917  case GEOS_MULTIPOLYGON:
918  {
919  QgsMultiPolygonV2* multiPolygon = new QgsMultiPolygonV2();
920 
921  int nParts = GEOSGetNumGeometries_r( geosinit.ctxt, geos );
922  for ( int i = 0; i < nParts; ++i )
923  {
924  QgsPolygonV2* poly = fromGeosPolygon( GEOSGetGeometryN_r( geosinit.ctxt, geos, i ) );
925  if ( poly )
926  {
927  multiPolygon->addGeometry( poly );
928  }
929  }
930  return multiPolygon;
931  }
932  case GEOS_GEOMETRYCOLLECTION:
933  {
934  QgsGeometryCollectionV2* geomCollection = new QgsGeometryCollectionV2();
935  int nParts = GEOSGetNumGeometries_r( geosinit.ctxt, geos );
936  for ( int i = 0; i < nParts; ++i )
937  {
938  QgsAbstractGeometryV2* geom = fromGeos( GEOSGetGeometryN_r( geosinit.ctxt, geos, i ) );
939  if ( geom )
940  {
941  geomCollection->addGeometry( geom );
942  }
943  }
944  return geomCollection;
945  }
946  }
947  return nullptr;
948 }
949 
950 QgsPolygonV2* QgsGeos::fromGeosPolygon( const GEOSGeometry* geos )
951 {
952  if ( GEOSGeomTypeId_r( geosinit.ctxt, geos ) != GEOS_POLYGON )
953  {
954  return nullptr;
955  }
956 
957  int nCoordDims = GEOSGeom_getCoordinateDimension_r( geosinit.ctxt, geos );
958  int nDims = GEOSGeom_getDimensions_r( geosinit.ctxt, geos );
959  bool hasZ = ( nCoordDims == 3 );
960  bool hasM = (( nDims - nCoordDims ) == 1 );
961 
962  QgsPolygonV2* polygon = new QgsPolygonV2();
963 
964  const GEOSGeometry* ring = GEOSGetExteriorRing_r( geosinit.ctxt, geos );
965  if ( ring )
966  {
967  polygon->setExteriorRing( sequenceToLinestring( ring, hasZ, hasM ) );
968  }
969 
970  QList<QgsCurveV2*> interiorRings;
971  for ( int i = 0; i < GEOSGetNumInteriorRings_r( geosinit.ctxt, geos ); ++i )
972  {
973  ring = GEOSGetInteriorRingN_r( geosinit.ctxt, geos, i );
974  if ( ring )
975  {
976  interiorRings.push_back( sequenceToLinestring( ring, hasZ, hasM ) );
977  }
978  }
979  polygon->setInteriorRings( interiorRings );
980 
981  return polygon;
982 }
983 
984 QgsLineStringV2* QgsGeos::sequenceToLinestring( const GEOSGeometry* geos, bool hasZ, bool hasM )
985 {
986  QgsPointSequenceV2 pts;
987  const GEOSCoordSequence* cs = GEOSGeom_getCoordSeq_r( geosinit.ctxt, geos );
988  unsigned int nPoints;
989  GEOSCoordSeq_getSize_r( geosinit.ctxt, cs, &nPoints );
990  pts.reserve( nPoints );
991  for ( unsigned int i = 0; i < nPoints; ++i )
992  {
993  pts.push_back( coordSeqPoint( cs, i, hasZ, hasM ) );
994  }
995  QgsLineStringV2* line = new QgsLineStringV2();
996  line->setPoints( pts );
997  return line;
998 }
999 
1000 int QgsGeos::numberOfGeometries( GEOSGeometry* g )
1001 {
1002  if ( !g )
1003  return 0;
1004 
1005  int geometryType = GEOSGeomTypeId_r( geosinit.ctxt, g );
1006  if ( geometryType == GEOS_POINT || geometryType == GEOS_LINESTRING || geometryType == GEOS_LINEARRING
1007  || geometryType == GEOS_POLYGON )
1008  return 1;
1009 
1010  //calling GEOSGetNumGeometries is save for multi types and collections also in geos2
1011  return GEOSGetNumGeometries_r( geosinit.ctxt, g );
1012 }
1013 
1014 QgsPointV2 QgsGeos::coordSeqPoint( const GEOSCoordSequence* cs, int i, bool hasZ, bool hasM )
1015 {
1016  if ( !cs )
1017  {
1018  return QgsPointV2();
1019  }
1020 
1021  double x, y;
1022  double z = 0;
1023  double m = 0;
1024  GEOSCoordSeq_getX_r( geosinit.ctxt, cs, i, &x );
1025  GEOSCoordSeq_getY_r( geosinit.ctxt, cs, i, &y );
1026  if ( hasZ )
1027  {
1028  GEOSCoordSeq_getZ_r( geosinit.ctxt, cs, i, &z );
1029  }
1030  if ( hasM )
1031  {
1032  GEOSCoordSeq_getOrdinate_r( geosinit.ctxt, cs, i, 3, &m );
1033  }
1034 
1036  if ( hasZ && hasM )
1037  {
1039  }
1040  else if ( hasZ )
1041  {
1042  t = QgsWKBTypes::PointZ;
1043  }
1044  else if ( hasM )
1045  {
1046  t = QgsWKBTypes::PointM;
1047  }
1048  return QgsPointV2( t, x, y, z, m );
1049 }
1050 
1051 GEOSGeometry* QgsGeos::asGeos( const QgsAbstractGeometryV2* geom, double precision )
1052 {
1053  int coordDims = 2;
1054  if ( geom->is3D() )
1055  {
1056  ++coordDims;
1057  }
1058  if ( geom->isMeasure() )
1059  {
1060  ++coordDims;
1061  }
1062 
1063  //curve
1064  const QgsCurveV2* curve = dynamic_cast< const QgsCurveV2* >( geom );
1065  const QgsCurvePolygonV2* curvePolygon = dynamic_cast< const QgsCurvePolygonV2* >( geom );
1066  if ( curve )
1067  {
1068  QScopedPointer< QgsLineStringV2> lineString( curve->curveToLine() );
1069  return createGeosLinestring( lineString.data(), precision );
1070  }
1071  else if ( curvePolygon )
1072  {
1073  QScopedPointer<QgsPolygonV2> polygon( curvePolygon->toPolygon() );
1074  return createGeosPolygon( polygon.data(), precision );
1075  }
1076  else if ( geom->geometryType() == "Point" )
1077  {
1078  return createGeosPoint( geom, coordDims, precision );
1079  }
1080  else if ( QgsWKBTypes::isMultiType( geom->wkbType() ) )
1081  {
1082  int geosType = GEOS_MULTIPOINT;
1083  if ( geom->geometryType() == "MultiPoint" )
1084  {
1085  geosType = GEOS_MULTIPOINT;
1086  }
1087  if ( geom->geometryType() == "MultiCurve" || geom->geometryType() == "MultiLineString" )
1088  {
1089  geosType = GEOS_MULTILINESTRING;
1090  }
1091  else if ( geom->geometryType() == "MultiSurface" || geom->geometryType() == "MultiPolygon" )
1092  {
1093  geosType = GEOS_MULTIPOLYGON;
1094  }
1095  else if ( geom->geometryType() == "GeometryCollection" )
1096  {
1097  geosType = GEOS_GEOMETRYCOLLECTION;
1098  }
1099 
1100  const QgsGeometryCollectionV2* c = dynamic_cast<const QgsGeometryCollectionV2*>( geom );
1101  if ( !c )
1102  {
1103  return nullptr;
1104  }
1105 
1106  QVector< GEOSGeometry* > geomVector( c->numGeometries() );
1107  for ( int i = 0; i < c->numGeometries(); ++i )
1108  {
1109  geomVector[i] = asGeos( c->geometryN( i ), precision );
1110  }
1111  return createGeosCollection( geosType, geomVector );
1112  }
1113 
1114  return nullptr;
1115 }
1116 
1117 QgsAbstractGeometryV2* QgsGeos::overlay( const QgsAbstractGeometryV2& geom, Overlay op, QString* errorMsg ) const
1118 {
1119  if ( !mGeos )
1120  {
1121  return nullptr;
1122  }
1123 
1124  GEOSGeomScopedPtr geosGeom( asGeos( &geom, mPrecision ) );
1125  if ( !geosGeom )
1126  {
1127  return nullptr;
1128  }
1129 
1130  try
1131  {
1132  GEOSGeomScopedPtr opGeom;
1133  switch ( op )
1134  {
1135  case INTERSECTION:
1136  opGeom.reset( GEOSIntersection_r( geosinit.ctxt, mGeos, geosGeom.get() ) );
1137  break;
1138  case DIFFERENCE:
1139  opGeom.reset( GEOSDifference_r( geosinit.ctxt, mGeos, geosGeom.get() ) );
1140  break;
1141  case UNION:
1142  {
1143  GEOSGeometry *unionGeometry = GEOSUnion_r( geosinit.ctxt, mGeos, geosGeom.get() );
1144 
1145  if ( unionGeometry && GEOSGeomTypeId_r( geosinit.ctxt, unionGeometry ) == GEOS_MULTILINESTRING )
1146  {
1147  GEOSGeometry *mergedLines = GEOSLineMerge_r( geosinit.ctxt, unionGeometry );
1148  if ( mergedLines )
1149  {
1150  GEOSGeom_destroy_r( geosinit.ctxt, unionGeometry );
1151  unionGeometry = mergedLines;
1152  }
1153  }
1154 
1155  opGeom.reset( unionGeometry );
1156  }
1157  break;
1158  case SYMDIFFERENCE:
1159  opGeom.reset( GEOSSymDifference_r( geosinit.ctxt, mGeos, geosGeom.get() ) );
1160  break;
1161  default: //unknown op
1162  return nullptr;
1163  }
1164  QgsAbstractGeometryV2* opResult = fromGeos( opGeom.get() );
1165  return opResult;
1166  }
1167  catch ( GEOSException &e )
1168  {
1169  if ( errorMsg )
1170  {
1171  *errorMsg = e.what();
1172  }
1173  return nullptr;
1174  }
1175 }
1176 
1177 bool QgsGeos::relation( const QgsAbstractGeometryV2& geom, Relation r, QString* errorMsg ) const
1178 {
1179  if ( !mGeos )
1180  {
1181  return false;
1182  }
1183 
1184  GEOSGeomScopedPtr geosGeom( asGeos( &geom, mPrecision ) );
1185  if ( !geosGeom )
1186  {
1187  return false;
1188  }
1189 
1190  bool result = false;
1191  try
1192  {
1193  if ( mGeosPrepared ) //use faster version with prepared geometry
1194  {
1195  switch ( r )
1196  {
1197  case INTERSECTS:
1198  result = ( GEOSPreparedIntersects_r( geosinit.ctxt, mGeosPrepared, geosGeom.get() ) == 1 );
1199  break;
1200  case TOUCHES:
1201  result = ( GEOSPreparedTouches_r( geosinit.ctxt, mGeosPrepared, geosGeom.get() ) == 1 );
1202  break;
1203  case CROSSES:
1204  result = ( GEOSPreparedCrosses_r( geosinit.ctxt, mGeosPrepared, geosGeom.get() ) == 1 );
1205  break;
1206  case WITHIN:
1207  result = ( GEOSPreparedWithin_r( geosinit.ctxt, mGeosPrepared, geosGeom.get() ) == 1 );
1208  break;
1209  case CONTAINS:
1210  result = ( GEOSPreparedContains_r( geosinit.ctxt, mGeosPrepared, geosGeom.get() ) == 1 );
1211  break;
1212  case DISJOINT:
1213  result = ( GEOSPreparedDisjoint_r( geosinit.ctxt, mGeosPrepared, geosGeom.get() ) == 1 );
1214  break;
1215  case OVERLAPS:
1216  result = ( GEOSPreparedOverlaps_r( geosinit.ctxt, mGeosPrepared, geosGeom.get() ) == 1 );
1217  break;
1218  default:
1219  return false;
1220  }
1221  return result;
1222  }
1223 
1224  switch ( r )
1225  {
1226  case INTERSECTS:
1227  result = ( GEOSIntersects_r( geosinit.ctxt, mGeos, geosGeom.get() ) == 1 );
1228  break;
1229  case TOUCHES:
1230  result = ( GEOSTouches_r( geosinit.ctxt, mGeos, geosGeom.get() ) == 1 );
1231  break;
1232  case CROSSES:
1233  result = ( GEOSCrosses_r( geosinit.ctxt, mGeos, geosGeom.get() ) == 1 );
1234  break;
1235  case WITHIN:
1236  result = ( GEOSWithin_r( geosinit.ctxt, mGeos, geosGeom.get() ) == 1 );
1237  break;
1238  case CONTAINS:
1239  result = ( GEOSContains_r( geosinit.ctxt, mGeos, geosGeom.get() ) == 1 );
1240  break;
1241  case DISJOINT:
1242  result = ( GEOSDisjoint_r( geosinit.ctxt, mGeos, geosGeom.get() ) == 1 );
1243  break;
1244  case OVERLAPS:
1245  result = ( GEOSOverlaps_r( geosinit.ctxt, mGeos, geosGeom.get() ) == 1 );
1246  break;
1247  default:
1248  return false;
1249  }
1250  }
1251  catch ( GEOSException &e )
1252  {
1253  if ( errorMsg )
1254  {
1255  *errorMsg = e.what();
1256  }
1257  return 0;
1258  }
1259 
1260  return result;
1261 }
1262 
1263 QgsAbstractGeometryV2* QgsGeos::buffer( double distance, int segments, QString* errorMsg ) const
1264 {
1265  if ( !mGeos )
1266  {
1267  return nullptr;
1268  }
1269 
1270  GEOSGeomScopedPtr geos;
1271  try
1272  {
1273  geos.reset( GEOSBuffer_r( geosinit.ctxt, mGeos, distance, segments ) );
1274  }
1275  CATCH_GEOS_WITH_ERRMSG( nullptr );
1276  return fromGeos( geos.get() );
1277 }
1278 
1279 QgsAbstractGeometryV2 *QgsGeos::buffer( double distance, int segments, int endCapStyle, int joinStyle, double mitreLimit, QString* errorMsg ) const
1280 {
1281  if ( !mGeos )
1282  {
1283  return nullptr;
1284  }
1285 
1286 #if defined(GEOS_VERSION_MAJOR) && defined(GEOS_VERSION_MINOR) && \
1287  ((GEOS_VERSION_MAJOR>3) || ((GEOS_VERSION_MAJOR==3) && (GEOS_VERSION_MINOR>=3)))
1288 
1289  GEOSGeomScopedPtr geos;
1290  try
1291  {
1292  geos.reset( GEOSBufferWithStyle_r( geosinit.ctxt, mGeos, distance, segments, endCapStyle, joinStyle, mitreLimit ) );
1293  }
1294  CATCH_GEOS_WITH_ERRMSG( nullptr );
1295  return fromGeos( geos.get() );
1296 #else
1297  return 0;
1298 #endif //0
1299 }
1300 
1301 QgsAbstractGeometryV2* QgsGeos::simplify( double tolerance, QString* errorMsg ) const
1302 {
1303  if ( !mGeos )
1304  {
1305  return nullptr;
1306  }
1307  GEOSGeomScopedPtr geos;
1308  try
1309  {
1310  geos.reset( GEOSTopologyPreserveSimplify_r( geosinit.ctxt, mGeos, tolerance ) );
1311  }
1312  CATCH_GEOS_WITH_ERRMSG( nullptr );
1313  return fromGeos( geos.get() );
1314 }
1315 
1317 {
1318  if ( !mGeos )
1319  {
1320  return nullptr;
1321  }
1322  GEOSGeomScopedPtr geos;
1323  try
1324  {
1325  geos.reset( GEOSInterpolate_r( geosinit.ctxt, mGeos, distance ) );
1326  }
1327  CATCH_GEOS_WITH_ERRMSG( nullptr );
1328  return fromGeos( geos.get() );
1329 }
1330 
1331 bool QgsGeos::centroid( QgsPointV2& pt, QString* errorMsg ) const
1332 {
1333  if ( !mGeos )
1334  {
1335  return false;
1336  }
1337 
1338  GEOSGeomScopedPtr geos;
1339  try
1340  {
1341  geos.reset( GEOSGetCentroid_r( geosinit.ctxt, mGeos ) );
1342  }
1343  CATCH_GEOS_WITH_ERRMSG( false );
1344 
1345  if ( !geos )
1346  {
1347  return false;
1348  }
1349 
1350  double x, y;
1351  GEOSGeomGetX_r( geosinit.ctxt, geos.get(), &x );
1352  GEOSGeomGetY_r( geosinit.ctxt, geos.get(), &y );
1353  pt.setX( x );
1354  pt.setY( y );
1355  return true;
1356 }
1357 
1359 {
1360  if ( !mGeos )
1361  {
1362  return nullptr;
1363  }
1364  GEOSGeomScopedPtr geos;
1365  try
1366  {
1367  geos.reset( GEOSEnvelope_r( geosinit.ctxt, mGeos ) );
1368  }
1369  CATCH_GEOS_WITH_ERRMSG( nullptr );
1370  return fromGeos( geos.get() );
1371 }
1372 
1373 bool QgsGeos::pointOnSurface( QgsPointV2& pt, QString* errorMsg ) const
1374 {
1375  if ( !mGeos )
1376  {
1377  return false;
1378  }
1379 
1380  GEOSGeomScopedPtr geos;
1381  try
1382  {
1383  geos.reset( GEOSPointOnSurface_r( geosinit.ctxt, mGeos ) );
1384 
1385  if ( !geos || GEOSisEmpty_r( geosinit.ctxt, geos.get() ) != 0 )
1386  {
1387  return false;
1388  }
1389 
1390  double x, y;
1391  GEOSGeomGetX_r( geosinit.ctxt, geos.get(), &x );
1392  GEOSGeomGetY_r( geosinit.ctxt, geos.get(), &y );
1393 
1394  pt.setX( x );
1395  pt.setY( y );
1396  }
1397  CATCH_GEOS_WITH_ERRMSG( false );
1398 
1399  return true;
1400 }
1401 
1403 {
1404  if ( !mGeos )
1405  {
1406  return nullptr;
1407  }
1408 
1409  try
1410  {
1411  GEOSGeometry* cHull = GEOSConvexHull_r( geosinit.ctxt, mGeos );
1412  QgsAbstractGeometryV2* cHullGeom = fromGeos( cHull );
1413  GEOSGeom_destroy_r( geosinit.ctxt, cHull );
1414  return cHullGeom;
1415  }
1416  CATCH_GEOS_WITH_ERRMSG( nullptr );
1417 }
1418 
1419 bool QgsGeos::isValid( QString* errorMsg ) const
1420 {
1421  if ( !mGeos )
1422  {
1423  return false;
1424  }
1425 
1426  try
1427  {
1428  return GEOSisValid_r( geosinit.ctxt, mGeos );
1429  }
1430  CATCH_GEOS_WITH_ERRMSG( false );
1431 }
1432 
1433 bool QgsGeos::isEqual( const QgsAbstractGeometryV2& geom, QString* errorMsg ) const
1434 {
1435  if ( !mGeos )
1436  {
1437  return false;
1438  }
1439 
1440  try
1441  {
1442  GEOSGeomScopedPtr geosGeom( asGeos( &geom, mPrecision ) );
1443  if ( !geosGeom )
1444  {
1445  return false;
1446  }
1447  bool equal = GEOSEquals_r( geosinit.ctxt, mGeos, geosGeom.get() );
1448  return equal;
1449  }
1450  CATCH_GEOS_WITH_ERRMSG( false );
1451 }
1452 
1453 bool QgsGeos::isEmpty( QString* errorMsg ) const
1454 {
1455  if ( !mGeos )
1456  {
1457  return false;
1458  }
1459 
1460  try
1461  {
1462  return GEOSisEmpty_r( geosinit.ctxt, mGeos );
1463  }
1464  CATCH_GEOS_WITH_ERRMSG( false );
1465 }
1466 
1467 GEOSCoordSequence* QgsGeos::createCoordinateSequence( const QgsCurveV2* curve, double precision )
1468 {
1469  bool segmentize = false;
1470  const QgsLineStringV2* line = dynamic_cast<const QgsLineStringV2*>( curve );
1471 
1472  if ( !line )
1473  {
1474  line = curve->curveToLine();
1475  segmentize = true;
1476  }
1477 
1478  if ( !line )
1479  {
1480  return nullptr;
1481  }
1482 
1483  bool hasZ = line->is3D();
1484  bool hasM = false; //line->isMeasure(); //disabled until geos supports m-coordinates
1485  int coordDims = 2;
1486  if ( hasZ )
1487  {
1488  ++coordDims;
1489  }
1490  if ( hasM )
1491  {
1492  ++coordDims;
1493  }
1494 
1495  int numPoints = line->numPoints();
1496  GEOSCoordSequence* coordSeq = nullptr;
1497  try
1498  {
1499  coordSeq = GEOSCoordSeq_create_r( geosinit.ctxt, numPoints, coordDims );
1500  if ( !coordSeq )
1501  {
1502  QgsMessageLog::logMessage( QObject::tr( "Could not create coordinate sequence for %1 points in %2 dimensions" ).arg( numPoints ).arg( coordDims ), QObject::tr( "GEOS" ) );
1503  return nullptr;
1504  }
1505  if ( precision > 0. )
1506  {
1507  for ( int i = 0; i < numPoints; ++i )
1508  {
1509  QgsPointV2 pt = line->pointN( i ); //todo: create method to get const point reference
1510  GEOSCoordSeq_setX_r( geosinit.ctxt, coordSeq, i, qgsRound( pt.x() / precision ) * precision );
1511  GEOSCoordSeq_setY_r( geosinit.ctxt, coordSeq, i, qgsRound( pt.y() / precision ) * precision );
1512  if ( hasZ )
1513  {
1514  GEOSCoordSeq_setOrdinate_r( geosinit.ctxt, coordSeq, i, 2, qgsRound( pt.z() / precision ) * precision );
1515  }
1516  if ( hasM )
1517  {
1518  GEOSCoordSeq_setOrdinate_r( geosinit.ctxt, coordSeq, i, 3, pt.m() );
1519  }
1520  }
1521  }
1522  else
1523  {
1524  for ( int i = 0; i < numPoints; ++i )
1525  {
1526  QgsPointV2 pt = line->pointN( i ); //todo: create method to get const point reference
1527  GEOSCoordSeq_setX_r( geosinit.ctxt, coordSeq, i, pt.x() );
1528  GEOSCoordSeq_setY_r( geosinit.ctxt, coordSeq, i, pt.y() );
1529  if ( hasZ )
1530  {
1531  GEOSCoordSeq_setOrdinate_r( geosinit.ctxt, coordSeq, i, 2, pt.z() );
1532  }
1533  if ( hasM )
1534  {
1535  GEOSCoordSeq_setOrdinate_r( geosinit.ctxt, coordSeq, i, 3, pt.m() );
1536  }
1537  }
1538  }
1539  }
1540  CATCH_GEOS( nullptr )
1541 
1542  if ( segmentize )
1543  {
1544  delete line;
1545  }
1546  return coordSeq;
1547 }
1548 
1549 GEOSGeometry* QgsGeos::createGeosPoint( const QgsAbstractGeometryV2* point, int coordDims, double precision )
1550 {
1551  const QgsPointV2* pt = dynamic_cast<const QgsPointV2*>( point );
1552  if ( !pt )
1553  return nullptr;
1554 
1555  GEOSGeometry* geosPoint = nullptr;
1556 
1557  try
1558  {
1559  GEOSCoordSequence* coordSeq = GEOSCoordSeq_create_r( geosinit.ctxt, 1, coordDims );
1560  if ( !coordSeq )
1561  {
1562  QgsMessageLog::logMessage( QObject::tr( "Could not create coordinate sequence for point with %1 dimensions" ).arg( coordDims ), QObject::tr( "GEOS" ) );
1563  return nullptr;
1564  }
1565  if ( precision > 0. )
1566  {
1567  GEOSCoordSeq_setX_r( geosinit.ctxt, coordSeq, 0, qgsRound( pt->x() / precision ) * precision );
1568  GEOSCoordSeq_setY_r( geosinit.ctxt, coordSeq, 0, qgsRound( pt->y() / precision ) * precision );
1569  if ( pt->is3D() )
1570  {
1571  GEOSCoordSeq_setOrdinate_r( geosinit.ctxt, coordSeq, 0, 2, qgsRound( pt->z() / precision ) * precision );
1572  }
1573  }
1574  else
1575  {
1576  GEOSCoordSeq_setX_r( geosinit.ctxt, coordSeq, 0, pt->x() );
1577  GEOSCoordSeq_setY_r( geosinit.ctxt, coordSeq, 0, pt->y() );
1578  if ( pt->is3D() )
1579  {
1580  GEOSCoordSeq_setOrdinate_r( geosinit.ctxt, coordSeq, 0, 2, pt->z() );
1581  }
1582  }
1583 #if 0 //disabled until geos supports m-coordinates
1584  if ( pt->isMeasure() )
1585  {
1586  GEOSCoordSeq_setOrdinate_r( geosinit.ctxt, coordSeq, 0, 3, pt->m() );
1587  }
1588 #endif
1589  geosPoint = GEOSGeom_createPoint_r( geosinit.ctxt, coordSeq );
1590  }
1591  CATCH_GEOS( nullptr )
1592  return geosPoint;
1593 }
1594 
1595 GEOSGeometry* QgsGeos::createGeosLinestring( const QgsAbstractGeometryV2* curve , double precision )
1596 {
1597  const QgsCurveV2* c = dynamic_cast<const QgsCurveV2*>( curve );
1598  if ( !c )
1599  return nullptr;
1600 
1601  GEOSCoordSequence* coordSeq = createCoordinateSequence( c, precision );
1602  if ( !coordSeq )
1603  return nullptr;
1604 
1605  GEOSGeometry* geosGeom = nullptr;
1606  try
1607  {
1608  geosGeom = GEOSGeom_createLineString_r( geosinit.ctxt, coordSeq );
1609  }
1610  CATCH_GEOS( nullptr )
1611  return geosGeom;
1612 }
1613 
1614 GEOSGeometry* QgsGeos::createGeosPolygon( const QgsAbstractGeometryV2* poly , double precision )
1615 {
1616  const QgsCurvePolygonV2* polygon = dynamic_cast<const QgsCurvePolygonV2*>( poly );
1617  if ( !polygon )
1618  return nullptr;
1619 
1620  const QgsCurveV2* exteriorRing = polygon->exteriorRing();
1621  if ( !exteriorRing )
1622  {
1623  return nullptr;
1624  }
1625 
1626  GEOSGeometry* geosPolygon = nullptr;
1627  try
1628  {
1629  GEOSGeometry* exteriorRingGeos = GEOSGeom_createLinearRing_r( geosinit.ctxt, createCoordinateSequence( exteriorRing, precision ) );
1630 
1631 
1632  int nHoles = polygon->numInteriorRings();
1633  GEOSGeometry** holes = nullptr;
1634  if ( nHoles > 0 )
1635  {
1636  holes = new GEOSGeometry*[ nHoles ];
1637  }
1638 
1639  for ( int i = 0; i < nHoles; ++i )
1640  {
1641  const QgsCurveV2* interiorRing = polygon->interiorRing( i );
1642  holes[i] = GEOSGeom_createLinearRing_r( geosinit.ctxt, createCoordinateSequence( interiorRing, precision ) );
1643  }
1644  geosPolygon = GEOSGeom_createPolygon_r( geosinit.ctxt, exteriorRingGeos, holes, nHoles );
1645  delete[] holes;
1646  }
1647  CATCH_GEOS( nullptr )
1648 
1649  return geosPolygon;
1650 }
1651 
1652 QgsAbstractGeometryV2* QgsGeos::offsetCurve( double distance, int segments, int joinStyle, double mitreLimit, QString* errorMsg ) const
1653 {
1654  if ( !mGeos )
1655  return nullptr;
1656 
1657  GEOSGeometry* offset = nullptr;
1658  try
1659  {
1660  offset = GEOSOffsetCurve_r( geosinit.ctxt, mGeos, distance, segments, joinStyle, mitreLimit );
1661  }
1662  CATCH_GEOS_WITH_ERRMSG( nullptr )
1663  QgsAbstractGeometryV2* offsetGeom = fromGeos( offset );
1664  GEOSGeom_destroy_r( geosinit.ctxt, offset );
1665  return offsetGeom;
1666 }
1667 
1668 QgsAbstractGeometryV2* QgsGeos::reshapeGeometry( const QgsLineStringV2& reshapeWithLine, int* errorCode, QString* errorMsg ) const
1669 {
1670  if ( !mGeos || reshapeWithLine.numPoints() < 2 || mGeometry->dimension() == 0 )
1671  {
1672  if ( errorCode ) { *errorCode = 1; }
1673  return nullptr;
1674  }
1675 
1676  GEOSGeometry* reshapeLineGeos = createGeosLinestring( &reshapeWithLine, mPrecision );
1677 
1678  //single or multi?
1679  int numGeoms = GEOSGetNumGeometries_r( geosinit.ctxt, mGeos );
1680  if ( numGeoms == -1 )
1681  {
1682  if ( errorCode ) { *errorCode = 1; }
1683  GEOSGeom_destroy_r( geosinit.ctxt, reshapeLineGeos );
1684  return nullptr;
1685  }
1686 
1687  bool isMultiGeom = false;
1688  int geosTypeId = GEOSGeomTypeId_r( geosinit.ctxt, mGeos );
1689  if ( geosTypeId == GEOS_MULTILINESTRING || geosTypeId == GEOS_MULTIPOLYGON )
1690  isMultiGeom = true;
1691 
1692  bool isLine = ( mGeometry->dimension() == 1 );
1693 
1694  if ( !isMultiGeom )
1695  {
1696  GEOSGeometry* reshapedGeometry;
1697  if ( isLine )
1698  {
1699  reshapedGeometry = reshapeLine( mGeos, reshapeLineGeos, mPrecision );
1700  }
1701  else
1702  {
1703  reshapedGeometry = reshapePolygon( mGeos, reshapeLineGeos, mPrecision );
1704  }
1705 
1706  if ( errorCode ) { *errorCode = 0; }
1707  QgsAbstractGeometryV2* reshapeResult = fromGeos( reshapedGeometry );
1708  GEOSGeom_destroy_r( geosinit.ctxt, reshapedGeometry );
1709  GEOSGeom_destroy_r( geosinit.ctxt, reshapeLineGeos );
1710  return reshapeResult;
1711  }
1712  else
1713  {
1714  try
1715  {
1716  //call reshape for each geometry part and replace mGeos with new geometry if reshape took place
1717  bool reshapeTookPlace = false;
1718 
1719  GEOSGeometry* currentReshapeGeometry = nullptr;
1720  GEOSGeometry** newGeoms = new GEOSGeometry*[numGeoms];
1721 
1722  for ( int i = 0; i < numGeoms; ++i )
1723  {
1724  if ( isLine )
1725  currentReshapeGeometry = reshapeLine( GEOSGetGeometryN_r( geosinit.ctxt, mGeos, i ), reshapeLineGeos, mPrecision );
1726  else
1727  currentReshapeGeometry = reshapePolygon( GEOSGetGeometryN_r( geosinit.ctxt, mGeos, i ), reshapeLineGeos, mPrecision );
1728 
1729  if ( currentReshapeGeometry )
1730  {
1731  newGeoms[i] = currentReshapeGeometry;
1732  reshapeTookPlace = true;
1733  }
1734  else
1735  {
1736  newGeoms[i] = GEOSGeom_clone_r( geosinit.ctxt, GEOSGetGeometryN_r( geosinit.ctxt, mGeos, i ) );
1737  }
1738  }
1739  GEOSGeom_destroy_r( geosinit.ctxt, reshapeLineGeos );
1740 
1741  GEOSGeometry* newMultiGeom = nullptr;
1742  if ( isLine )
1743  {
1744  newMultiGeom = GEOSGeom_createCollection_r( geosinit.ctxt, GEOS_MULTILINESTRING, newGeoms, numGeoms );
1745  }
1746  else //multipolygon
1747  {
1748  newMultiGeom = GEOSGeom_createCollection_r( geosinit.ctxt, GEOS_MULTIPOLYGON, newGeoms, numGeoms );
1749  }
1750 
1751  delete[] newGeoms;
1752  if ( !newMultiGeom )
1753  {
1754  if ( errorCode ) { *errorCode = 3; }
1755  return nullptr;
1756  }
1757 
1758  if ( reshapeTookPlace )
1759  {
1760  if ( errorCode ) { *errorCode = 0; }
1761  QgsAbstractGeometryV2* reshapedMultiGeom = fromGeos( newMultiGeom );
1762  GEOSGeom_destroy_r( geosinit.ctxt, newMultiGeom );
1763  return reshapedMultiGeom;
1764  }
1765  else
1766  {
1767  GEOSGeom_destroy_r( geosinit.ctxt, newMultiGeom );
1768  if ( errorCode ) { *errorCode = 1; }
1769  return nullptr;
1770  }
1771  }
1772  CATCH_GEOS_WITH_ERRMSG( nullptr )
1773  }
1774 }
1775 
1776 QgsGeometry QgsGeos::closestPoint( const QgsGeometry& other, QString* errorMsg ) const
1777 {
1778  if ( !mGeos || other.isEmpty() )
1779  {
1780  return QgsGeometry();
1781  }
1782 
1783  GEOSGeomScopedPtr otherGeom( asGeos( other.geometry(), mPrecision ) );
1784  if ( !otherGeom )
1785  {
1786  return QgsGeometry();
1787  }
1788 
1789  double nx = 0.0;
1790  double ny = 0.0;
1791  try
1792  {
1793  GEOSCoordSequence* nearestCoord = GEOSNearestPoints_r( geosinit.ctxt, mGeos, otherGeom.get() );
1794 
1795  ( void )GEOSCoordSeq_getX_r( geosinit.ctxt, nearestCoord, 0, &nx );
1796  ( void )GEOSCoordSeq_getY_r( geosinit.ctxt, nearestCoord, 0, &ny );
1797  GEOSCoordSeq_destroy_r( geosinit.ctxt, nearestCoord );
1798  }
1799  catch ( GEOSException &e )
1800  {
1801  if ( errorMsg )
1802  {
1803  *errorMsg = e.what();
1804  }
1805  return QgsGeometry();
1806  }
1807 
1808  return QgsGeometry( new QgsPointV2( nx, ny ) );
1809 }
1810 
1811 QgsGeometry QgsGeos::shortestLine( const QgsGeometry& other, QString* errorMsg ) const
1812 {
1813  if ( !mGeos || other.isEmpty() )
1814  {
1815  return QgsGeometry();
1816  }
1817 
1818  GEOSGeomScopedPtr otherGeom( asGeos( other.geometry(), mPrecision ) );
1819  if ( !otherGeom )
1820  {
1821  return QgsGeometry();
1822  }
1823 
1824  double nx1 = 0.0;
1825  double ny1 = 0.0;
1826  double nx2 = 0.0;
1827  double ny2 = 0.0;
1828  try
1829  {
1830  GEOSCoordSequence* nearestCoord = GEOSNearestPoints_r( geosinit.ctxt, mGeos, otherGeom.get() );
1831 
1832  ( void )GEOSCoordSeq_getX_r( geosinit.ctxt, nearestCoord, 0, &nx1 );
1833  ( void )GEOSCoordSeq_getY_r( geosinit.ctxt, nearestCoord, 0, &ny1 );
1834  ( void )GEOSCoordSeq_getX_r( geosinit.ctxt, nearestCoord, 1, &nx2 );
1835  ( void )GEOSCoordSeq_getY_r( geosinit.ctxt, nearestCoord, 1, &ny2 );
1836 
1837  GEOSCoordSeq_destroy_r( geosinit.ctxt, nearestCoord );
1838  }
1839  catch ( GEOSException &e )
1840  {
1841  if ( errorMsg )
1842  {
1843  *errorMsg = e.what();
1844  }
1845  return QgsGeometry();
1846  }
1847 
1848  QgsLineStringV2* line = new QgsLineStringV2();
1849  line->addVertex( QgsPointV2( nx1, ny1 ) );
1850  line->addVertex( QgsPointV2( nx2, ny2 ) );
1851  return QgsGeometry( line );
1852 }
1853 
1854 GEOSGeometry* QgsGeos::reshapeLine( const GEOSGeometry* line, const GEOSGeometry* reshapeLineGeos , double precision )
1855 {
1856  if ( !line || !reshapeLineGeos )
1857  return nullptr;
1858 
1859  bool atLeastTwoIntersections = false;
1860 
1861  try
1862  {
1863  //make sure there are at least two intersection between line and reshape geometry
1864  GEOSGeometry* intersectGeom = GEOSIntersection_r( geosinit.ctxt, line, reshapeLineGeos );
1865  if ( intersectGeom )
1866  {
1867  atLeastTwoIntersections = ( GEOSGeomTypeId_r( geosinit.ctxt, intersectGeom ) == GEOS_MULTIPOINT
1868  && GEOSGetNumGeometries_r( geosinit.ctxt, intersectGeom ) > 1 );
1869  GEOSGeom_destroy_r( geosinit.ctxt, intersectGeom );
1870  }
1871  }
1872  catch ( GEOSException &e )
1873  {
1874  QgsMessageLog::logMessage( QObject::tr( "Exception: %1" ).arg( e.what() ), QObject::tr( "GEOS" ) );
1875  atLeastTwoIntersections = false;
1876  }
1877 
1878  if ( !atLeastTwoIntersections )
1879  return nullptr;
1880 
1881  //begin and end point of original line
1882  const GEOSCoordSequence* lineCoordSeq = GEOSGeom_getCoordSeq_r( geosinit.ctxt, line );
1883  if ( !lineCoordSeq )
1884  return nullptr;
1885 
1886  unsigned int lineCoordSeqSize;
1887  if ( GEOSCoordSeq_getSize_r( geosinit.ctxt, lineCoordSeq, &lineCoordSeqSize ) == 0 )
1888  return nullptr;
1889 
1890  if ( lineCoordSeqSize < 2 )
1891  return nullptr;
1892 
1893  //first and last vertex of line
1894  double x1, y1, x2, y2;
1895  GEOSCoordSeq_getX_r( geosinit.ctxt, lineCoordSeq, 0, &x1 );
1896  GEOSCoordSeq_getY_r( geosinit.ctxt, lineCoordSeq, 0, &y1 );
1897  GEOSCoordSeq_getX_r( geosinit.ctxt, lineCoordSeq, lineCoordSeqSize - 1, &x2 );
1898  GEOSCoordSeq_getY_r( geosinit.ctxt, lineCoordSeq, lineCoordSeqSize - 1, &y2 );
1899  QgsPointV2 beginPoint( x1, y1 );
1900  GEOSGeometry* beginLineVertex = createGeosPoint( &beginPoint, 2, precision );
1901  QgsPointV2 endPoint( x2, y2 );
1902  GEOSGeometry* endLineVertex = createGeosPoint( &endPoint, 2, precision );
1903 
1904  bool isRing = false;
1905  if ( GEOSGeomTypeId_r( geosinit.ctxt, line ) == GEOS_LINEARRING
1906  || GEOSEquals_r( geosinit.ctxt, beginLineVertex, endLineVertex ) == 1 )
1907  isRing = true;
1908 
1909  //node line and reshape line
1910  GEOSGeometry* nodedGeometry = nodeGeometries( reshapeLineGeos, line );
1911  if ( !nodedGeometry )
1912  {
1913  GEOSGeom_destroy_r( geosinit.ctxt, beginLineVertex );
1914  GEOSGeom_destroy_r( geosinit.ctxt, endLineVertex );
1915  return nullptr;
1916  }
1917 
1918  //and merge them together
1919  GEOSGeometry *mergedLines = GEOSLineMerge_r( geosinit.ctxt, nodedGeometry );
1920  GEOSGeom_destroy_r( geosinit.ctxt, nodedGeometry );
1921  if ( !mergedLines )
1922  {
1923  GEOSGeom_destroy_r( geosinit.ctxt, beginLineVertex );
1924  GEOSGeom_destroy_r( geosinit.ctxt, endLineVertex );
1925  return nullptr;
1926  }
1927 
1928  int numMergedLines = GEOSGetNumGeometries_r( geosinit.ctxt, mergedLines );
1929  if ( numMergedLines < 2 ) //some special cases. Normally it is >2
1930  {
1931  GEOSGeom_destroy_r( geosinit.ctxt, beginLineVertex );
1932  GEOSGeom_destroy_r( geosinit.ctxt, endLineVertex );
1933  if ( numMergedLines == 1 ) //reshape line is from begin to endpoint. So we keep the reshapeline
1934  return GEOSGeom_clone_r( geosinit.ctxt, reshapeLineGeos );
1935  else
1936  return nullptr;
1937  }
1938 
1939  QList<GEOSGeometry*> resultLineParts; //collection with the line segments that will be contained in result
1940  QList<GEOSGeometry*> probableParts; //parts where we can decide on inclusion only after going through all the candidates
1941 
1942  for ( int i = 0; i < numMergedLines; ++i )
1943  {
1944  const GEOSGeometry* currentGeom;
1945 
1946  currentGeom = GEOSGetGeometryN_r( geosinit.ctxt, mergedLines, i );
1947  const GEOSCoordSequence* currentCoordSeq = GEOSGeom_getCoordSeq_r( geosinit.ctxt, currentGeom );
1948  unsigned int currentCoordSeqSize;
1949  GEOSCoordSeq_getSize_r( geosinit.ctxt, currentCoordSeq, &currentCoordSeqSize );
1950  if ( currentCoordSeqSize < 2 )
1951  continue;
1952 
1953  //get the two endpoints of the current line merge result
1954  double xBegin, xEnd, yBegin, yEnd;
1955  GEOSCoordSeq_getX_r( geosinit.ctxt, currentCoordSeq, 0, &xBegin );
1956  GEOSCoordSeq_getY_r( geosinit.ctxt, currentCoordSeq, 0, &yBegin );
1957  GEOSCoordSeq_getX_r( geosinit.ctxt, currentCoordSeq, currentCoordSeqSize - 1, &xEnd );
1958  GEOSCoordSeq_getY_r( geosinit.ctxt, currentCoordSeq, currentCoordSeqSize - 1, &yEnd );
1959  QgsPointV2 beginPoint( xBegin, yBegin );
1960  GEOSGeometry* beginCurrentGeomVertex = createGeosPoint( &beginPoint, 2, precision );
1961  QgsPointV2 endPoint( xEnd, yEnd );
1962  GEOSGeometry* endCurrentGeomVertex = createGeosPoint( &endPoint, 2, precision );
1963 
1964  //check how many endpoints of the line merge result are on the (original) line
1965  int nEndpointsOnOriginalLine = 0;
1966  if ( pointContainedInLine( beginCurrentGeomVertex, line ) == 1 )
1967  nEndpointsOnOriginalLine += 1;
1968 
1969  if ( pointContainedInLine( endCurrentGeomVertex, line ) == 1 )
1970  nEndpointsOnOriginalLine += 1;
1971 
1972  //check how many endpoints equal the endpoints of the original line
1973  int nEndpointsSameAsOriginalLine = 0;
1974  if ( GEOSEquals_r( geosinit.ctxt, beginCurrentGeomVertex, beginLineVertex ) == 1
1975  || GEOSEquals_r( geosinit.ctxt, beginCurrentGeomVertex, endLineVertex ) == 1 )
1976  nEndpointsSameAsOriginalLine += 1;
1977 
1978  if ( GEOSEquals_r( geosinit.ctxt, endCurrentGeomVertex, beginLineVertex ) == 1
1979  || GEOSEquals_r( geosinit.ctxt, endCurrentGeomVertex, endLineVertex ) == 1 )
1980  nEndpointsSameAsOriginalLine += 1;
1981 
1982  //check if the current geometry overlaps the original geometry (GEOSOverlap does not seem to work with linestrings)
1983  bool currentGeomOverlapsOriginalGeom = false;
1984  bool currentGeomOverlapsReshapeLine = false;
1985  if ( lineContainedInLine( currentGeom, line ) == 1 )
1986  currentGeomOverlapsOriginalGeom = true;
1987 
1988  if ( lineContainedInLine( currentGeom, reshapeLineGeos ) == 1 )
1989  currentGeomOverlapsReshapeLine = true;
1990 
1991  //logic to decide if this part belongs to the result
1992  if ( nEndpointsSameAsOriginalLine == 1 && nEndpointsOnOriginalLine == 2 && currentGeomOverlapsOriginalGeom )
1993  {
1994  resultLineParts.push_back( GEOSGeom_clone_r( geosinit.ctxt, currentGeom ) );
1995  }
1996  //for closed rings, we take one segment from the candidate list
1997  else if ( isRing && nEndpointsOnOriginalLine == 2 && currentGeomOverlapsOriginalGeom )
1998  {
1999  probableParts.push_back( GEOSGeom_clone_r( geosinit.ctxt, currentGeom ) );
2000  }
2001  else if ( nEndpointsOnOriginalLine == 2 && !currentGeomOverlapsOriginalGeom )
2002  {
2003  resultLineParts.push_back( GEOSGeom_clone_r( geosinit.ctxt, currentGeom ) );
2004  }
2005  else if ( nEndpointsSameAsOriginalLine == 2 && !currentGeomOverlapsOriginalGeom )
2006  {
2007  resultLineParts.push_back( GEOSGeom_clone_r( geosinit.ctxt, currentGeom ) );
2008  }
2009  else if ( currentGeomOverlapsOriginalGeom && currentGeomOverlapsReshapeLine )
2010  {
2011  resultLineParts.push_back( GEOSGeom_clone_r( geosinit.ctxt, currentGeom ) );
2012  }
2013 
2014  GEOSGeom_destroy_r( geosinit.ctxt, beginCurrentGeomVertex );
2015  GEOSGeom_destroy_r( geosinit.ctxt, endCurrentGeomVertex );
2016  }
2017 
2018  //add the longest segment from the probable list for rings (only used for polygon rings)
2019  if ( isRing && !probableParts.isEmpty() )
2020  {
2021  GEOSGeometry* maxGeom = nullptr; //the longest geometry in the probabla list
2022  GEOSGeometry* currentGeom = nullptr;
2023  double maxLength = -std::numeric_limits<double>::max();
2024  double currentLength = 0;
2025  for ( int i = 0; i < probableParts.size(); ++i )
2026  {
2027  currentGeom = probableParts.at( i );
2028  GEOSLength_r( geosinit.ctxt, currentGeom, &currentLength );
2029  if ( currentLength > maxLength )
2030  {
2031  maxLength = currentLength;
2032  GEOSGeom_destroy_r( geosinit.ctxt, maxGeom );
2033  maxGeom = currentGeom;
2034  }
2035  else
2036  {
2037  GEOSGeom_destroy_r( geosinit.ctxt, currentGeom );
2038  }
2039  }
2040  resultLineParts.push_back( maxGeom );
2041  }
2042 
2043  GEOSGeom_destroy_r( geosinit.ctxt, beginLineVertex );
2044  GEOSGeom_destroy_r( geosinit.ctxt, endLineVertex );
2045  GEOSGeom_destroy_r( geosinit.ctxt, mergedLines );
2046 
2047  GEOSGeometry* result = nullptr;
2048  if ( resultLineParts.size() < 1 )
2049  return nullptr;
2050 
2051  if ( resultLineParts.size() == 1 ) //the whole result was reshaped
2052  {
2053  result = resultLineParts[0];
2054  }
2055  else //>1
2056  {
2057  GEOSGeometry **lineArray = new GEOSGeometry*[resultLineParts.size()];
2058  for ( int i = 0; i < resultLineParts.size(); ++i )
2059  {
2060  lineArray[i] = resultLineParts[i];
2061  }
2062 
2063  //create multiline from resultLineParts
2064  GEOSGeometry* multiLineGeom = GEOSGeom_createCollection_r( geosinit.ctxt, GEOS_MULTILINESTRING, lineArray, resultLineParts.size() );
2065  delete [] lineArray;
2066 
2067  //then do a linemerge with the newly combined partstrings
2068  result = GEOSLineMerge_r( geosinit.ctxt, multiLineGeom );
2069  GEOSGeom_destroy_r( geosinit.ctxt, multiLineGeom );
2070  }
2071 
2072  //now test if the result is a linestring. Otherwise something went wrong
2073  if ( GEOSGeomTypeId_r( geosinit.ctxt, result ) != GEOS_LINESTRING )
2074  {
2075  GEOSGeom_destroy_r( geosinit.ctxt, result );
2076  return nullptr;
2077  }
2078 
2079  return result;
2080 }
2081 
2082 GEOSGeometry* QgsGeos::reshapePolygon( const GEOSGeometry* polygon, const GEOSGeometry* reshapeLineGeos, double precision )
2083 {
2084  //go through outer shell and all inner rings and check if there is exactly one intersection of a ring and the reshape line
2085  int nIntersections = 0;
2086  int lastIntersectingRing = -2;
2087  const GEOSGeometry* lastIntersectingGeom = nullptr;
2088 
2089  int nRings = GEOSGetNumInteriorRings_r( geosinit.ctxt, polygon );
2090  if ( nRings < 0 )
2091  return nullptr;
2092 
2093  //does outer ring intersect?
2094  const GEOSGeometry* outerRing = GEOSGetExteriorRing_r( geosinit.ctxt, polygon );
2095  if ( GEOSIntersects_r( geosinit.ctxt, outerRing, reshapeLineGeos ) == 1 )
2096  {
2097  ++nIntersections;
2098  lastIntersectingRing = -1;
2099  lastIntersectingGeom = outerRing;
2100  }
2101 
2102  //do inner rings intersect?
2103  const GEOSGeometry **innerRings = new const GEOSGeometry*[nRings];
2104 
2105  try
2106  {
2107  for ( int i = 0; i < nRings; ++i )
2108  {
2109  innerRings[i] = GEOSGetInteriorRingN_r( geosinit.ctxt, polygon, i );
2110  if ( GEOSIntersects_r( geosinit.ctxt, innerRings[i], reshapeLineGeos ) == 1 )
2111  {
2112  ++nIntersections;
2113  lastIntersectingRing = i;
2114  lastIntersectingGeom = innerRings[i];
2115  }
2116  }
2117  }
2118  catch ( GEOSException &e )
2119  {
2120  QgsMessageLog::logMessage( QObject::tr( "Exception: %1" ).arg( e.what() ), QObject::tr( "GEOS" ) );
2121  nIntersections = 0;
2122  }
2123 
2124  if ( nIntersections != 1 ) //reshape line is only allowed to intersect one ring
2125  {
2126  delete [] innerRings;
2127  return nullptr;
2128  }
2129 
2130  //we have one intersecting ring, let's try to reshape it
2131  GEOSGeometry* reshapeResult = reshapeLine( lastIntersectingGeom, reshapeLineGeos, precision );
2132  if ( !reshapeResult )
2133  {
2134  delete [] innerRings;
2135  return nullptr;
2136  }
2137 
2138  //if reshaping took place, we need to reassemble the polygon and its rings
2139  GEOSGeometry* newRing = nullptr;
2140  const GEOSCoordSequence* reshapeSequence = GEOSGeom_getCoordSeq_r( geosinit.ctxt, reshapeResult );
2141  GEOSCoordSequence* newCoordSequence = GEOSCoordSeq_clone_r( geosinit.ctxt, reshapeSequence );
2142 
2143  GEOSGeom_destroy_r( geosinit.ctxt, reshapeResult );
2144 
2145  newRing = GEOSGeom_createLinearRing_r( geosinit.ctxt, newCoordSequence );
2146  if ( !newRing )
2147  {
2148  delete [] innerRings;
2149  return nullptr;
2150  }
2151 
2152  GEOSGeometry* newOuterRing = nullptr;
2153  if ( lastIntersectingRing == -1 )
2154  newOuterRing = newRing;
2155  else
2156  newOuterRing = GEOSGeom_clone_r( geosinit.ctxt, outerRing );
2157 
2158  //check if all the rings are still inside the outer boundary
2159  QList<GEOSGeometry*> ringList;
2160  if ( nRings > 0 )
2161  {
2162  GEOSGeometry* outerRingPoly = GEOSGeom_createPolygon_r( geosinit.ctxt, GEOSGeom_clone_r( geosinit.ctxt, newOuterRing ), nullptr, 0 );
2163  if ( outerRingPoly )
2164  {
2165  GEOSGeometry* currentRing = nullptr;
2166  for ( int i = 0; i < nRings; ++i )
2167  {
2168  if ( lastIntersectingRing == i )
2169  currentRing = newRing;
2170  else
2171  currentRing = GEOSGeom_clone_r( geosinit.ctxt, innerRings[i] );
2172 
2173  //possibly a ring is no longer contained in the result polygon after reshape
2174  if ( GEOSContains_r( geosinit.ctxt, outerRingPoly, currentRing ) == 1 )
2175  ringList.push_back( currentRing );
2176  else
2177  GEOSGeom_destroy_r( geosinit.ctxt, currentRing );
2178  }
2179  }
2180  GEOSGeom_destroy_r( geosinit.ctxt, outerRingPoly );
2181  }
2182 
2183  GEOSGeometry** newInnerRings = new GEOSGeometry*[ringList.size()];
2184  for ( int i = 0; i < ringList.size(); ++i )
2185  newInnerRings[i] = ringList.at( i );
2186 
2187  delete [] innerRings;
2188 
2189  GEOSGeometry* reshapedPolygon = GEOSGeom_createPolygon_r( geosinit.ctxt, newOuterRing, newInnerRings, ringList.size() );
2190  delete[] newInnerRings;
2191 
2192  return reshapedPolygon;
2193 }
2194 
2195 int QgsGeos::lineContainedInLine( const GEOSGeometry* line1, const GEOSGeometry* line2 )
2196 {
2197  if ( !line1 || !line2 )
2198  {
2199  return -1;
2200  }
2201 
2202  double bufferDistance = pow( 10.0L, geomDigits( line2 ) - 11 );
2203 
2204  GEOSGeometry* bufferGeom = GEOSBuffer_r( geosinit.ctxt, line2, bufferDistance, DEFAULT_QUADRANT_SEGMENTS );
2205  if ( !bufferGeom )
2206  return -2;
2207 
2208  GEOSGeometry* intersectionGeom = GEOSIntersection_r( geosinit.ctxt, bufferGeom, line1 );
2209 
2210  //compare ratio between line1Length and intersectGeomLength (usually close to 1 if line1 is contained in line2)
2211  double intersectGeomLength;
2212  double line1Length;
2213 
2214  GEOSLength_r( geosinit.ctxt, intersectionGeom, &intersectGeomLength );
2215  GEOSLength_r( geosinit.ctxt, line1, &line1Length );
2216 
2217  GEOSGeom_destroy_r( geosinit.ctxt, bufferGeom );
2218  GEOSGeom_destroy_r( geosinit.ctxt, intersectionGeom );
2219 
2220  double intersectRatio = line1Length / intersectGeomLength;
2221  if ( intersectRatio > 0.9 && intersectRatio < 1.1 )
2222  return 1;
2223 
2224  return 0;
2225 }
2226 
2227 int QgsGeos::pointContainedInLine( const GEOSGeometry* point, const GEOSGeometry* line )
2228 {
2229  if ( !point || !line )
2230  return -1;
2231 
2232  double bufferDistance = pow( 10.0L, geomDigits( line ) - 11 );
2233 
2234  GEOSGeometry* lineBuffer = GEOSBuffer_r( geosinit.ctxt, line, bufferDistance, 8 );
2235  if ( !lineBuffer )
2236  return -2;
2237 
2238  bool contained = false;
2239  if ( GEOSContains_r( geosinit.ctxt, lineBuffer, point ) == 1 )
2240  contained = true;
2241 
2242  GEOSGeom_destroy_r( geosinit.ctxt, lineBuffer );
2243  return contained;
2244 }
2245 
2246 int QgsGeos::geomDigits( const GEOSGeometry* geom )
2247 {
2248  GEOSGeomScopedPtr bbox( GEOSEnvelope_r( geosinit.ctxt, geom ) );
2249  if ( !bbox.get() )
2250  return -1;
2251 
2252  const GEOSGeometry* bBoxRing = GEOSGetExteriorRing_r( geosinit.ctxt, bbox.get() );
2253  if ( !bBoxRing )
2254  return -1;
2255 
2256  const GEOSCoordSequence* bBoxCoordSeq = GEOSGeom_getCoordSeq_r( geosinit.ctxt, bBoxRing );
2257 
2258  if ( !bBoxCoordSeq )
2259  return -1;
2260 
2261  unsigned int nCoords = 0;
2262  if ( !GEOSCoordSeq_getSize_r( geosinit.ctxt, bBoxCoordSeq, &nCoords ) )
2263  return -1;
2264 
2265  int maxDigits = -1;
2266  for ( unsigned int i = 0; i < nCoords - 1; ++i )
2267  {
2268  double t;
2269  GEOSCoordSeq_getX_r( geosinit.ctxt, bBoxCoordSeq, i, &t );
2270 
2271  int digits;
2272  digits = ceil( log10( fabs( t ) ) );
2273  if ( digits > maxDigits )
2274  maxDigits = digits;
2275 
2276  GEOSCoordSeq_getY_r( geosinit.ctxt, bBoxCoordSeq, i, &t );
2277  digits = ceil( log10( fabs( t ) ) );
2278  if ( digits > maxDigits )
2279  maxDigits = digits;
2280  }
2281 
2282  return maxDigits;
2283 }
2284 
2285 GEOSContextHandle_t QgsGeos::getGEOSHandler()
2286 {
2287  return geosinit.ctxt;
2288 }
virtual bool addGeometry(QgsAbstractGeometryV2 *g) override
Adds a geometry and takes ownership.
const QgsCurveV2 * exteriorRing() const
int splitGeometry(const QgsLineStringV2 &splitLine, QList< QgsAbstractGeometryV2 * > &newGeometries, bool topological, QgsPointSequenceV2 &topologyTestPoints, QString *errorMsg=nullptr) const override
Splits this geometry according to a given line.
Definition: qgsgeos.cpp:378
void clear()
void reset(GEOSGeometry *geom)
Definition: qgsgeos.cpp:121
QgsWKBTypes::Type wkbType() const
Returns the WKB type of the geometry.
QgsPointV2 pointN(int i) const
Returns the specified point from inside the line string.
#define CATCH_GEOS(r)
Definition: qgsgeos.cpp:34
QgsAbstractGeometryV2 * combine(const QgsAbstractGeometryV2 &geom, QString *errorMsg=nullptr) const override
Definition: qgsgeos.cpp:188
static QgsPointV2 coordSeqPoint(const GEOSCoordSequence *cs, int i, bool hasZ, bool hasM)
Definition: qgsgeos.cpp:1014
double x() const
Returns the point&#39;s x-coordinate.
Definition: qgspointv2.h:68
void push_back(const T &value)
bool isEmpty(QString *errorMsg=nullptr) const override
Definition: qgsgeos.cpp:1453
virtual bool addGeometry(QgsAbstractGeometryV2 *g) override
Adds a geometry and takes ownership.
#define QgsDebugMsg(str)
Definition: qgslogger.h:33
bool relatePattern(const QgsAbstractGeometryV2 &geom, const QString &pattern, QString *errorMsg=nullptr) const override
Tests whether two geometries are related by a specified Dimensional Extended 9 Intersection Model (DE...
Definition: qgsgeos.cpp:316
QgsAbstractGeometryV2 * geometry() const
Returns the underlying geometry store.
Multi curve geometry collection.
void reserve(int alloc)
double area(QString *errorMsg=nullptr) const override
Definition: qgsgeos.cpp:345
const_iterator constEnd() const
const T & at(int i) const
static bool isMultiType(Type type)
Returns true if the WKB type is a multi type.
Definition: qgswkbtypes.h:476
QgsGeometry closestPoint(const QgsGeometry &other, QString *errorMsg=nullptr) const
Returns the closest point on the geometry to the other geometry.
Definition: qgsgeos.cpp:1776
Abstract base class for all geometries.
double length(QString *errorMsg=nullptr) const override
Definition: qgsgeos.cpp:362
A geometry is the spatial representation of a feature.
Definition: qgsgeometry.h:76
static QgsAbstractGeometryV2 * fromGeos(const GEOSGeometry *geos)
Create a geometry from a GEOSGeometry.
Definition: qgsgeos.cpp:862
void setX(double x)
Sets the point&#39;s x-coordinate.
Definition: qgspointv2.h:124
QgsAbstractGeometryV2 * symDifference(const QgsAbstractGeometryV2 &geom, QString *errorMsg=nullptr) const override
Definition: qgsgeos.cpp:217
Scoped GEOS pointer.
Definition: qgsgeos.cpp:114
~QgsGeos()
Definition: qgsgeos.cpp:141
bool isValid(QString *errorMsg=nullptr) const override
Definition: qgsgeos.cpp:1419
QgsGeos(const QgsAbstractGeometryV2 *geometry, double precision=0)
GEOS geometry engine constructor.
Definition: qgsgeos.cpp:135
Multi point geometry collection.
QgsAbstractGeometryV2 * simplify(double tolerance, QString *errorMsg=nullptr) const override
Definition: qgsgeos.cpp:1301
virtual bool addGeometry(QgsAbstractGeometryV2 *g) override
Adds a geometry and takes ownership.
Multi line string geometry collection.
QString tr(const char *sourceText, const char *disambiguation, int n)
void setY(double y)
Sets the point&#39;s y-coordinate.
Definition: qgspointv2.h:130
static GEOSContextHandle_t getGEOSHandler()
Definition: qgsgeos.cpp:2285
int size() const
GEOSGeometry * get() const
Definition: qgsgeos.cpp:119
QgsAbstractGeometryV2 * interpolate(double distance, QString *errorMsg=nullptr) const override
Definition: qgsgeos.cpp:1316
GEOSGeomScopedPtr(GEOSGeometry *geom=nullptr)
Definition: qgsgeos.cpp:117
double y() const
Returns the point&#39;s y-coordinate.
Definition: qgspointv2.h:74
bool within(const QgsAbstractGeometryV2 &geom, QString *errorMsg=nullptr) const override
Definition: qgsgeos.cpp:262
bool contains(const QgsAbstractGeometryV2 &geom, QString *errorMsg=nullptr) const override
Definition: qgsgeos.cpp:272
double ANALYSIS_EXPORT max(double x, double y)
Returns the maximum of two doubles or the first argument if both are equal.
bool intersects(const QgsAbstractGeometryV2 &geom, QString *errorMsg=nullptr) const override
Definition: qgsgeos.cpp:247
Polygon geometry type.
Definition: qgspolygonv2.h:29
#define CATCH_GEOS_WITH_ERRMSG(r)
Definition: qgsgeos.cpp:41
double qgsRound(double x)
Definition: qgis.h:312
void clear()
bool is3D() const
Returns true if the geometry is 3D and contains a z-value.
void setInteriorRings(const QList< QgsCurveV2 * > &rings)
Sets all interior rings (takes ownership)
QString fromUtf8(const char *str, int size)
virtual void setExteriorRing(QgsCurveV2 *ring) override
Sets the exterior ring of the polygon.
void resize(int size)
static QgsPolygonV2 * fromGeosPolygon(const GEOSGeometry *geos)
Definition: qgsgeos.cpp:950
bool isMeasure() const
Returns true if the geometry contains m values.
double z() const
Returns the point&#39;s z-coordinate.
Definition: qgspointv2.h:80
Line string geometry type, with support for z-dimension and m-values.
void setPoints(const QgsPointSequenceV2 &points)
Resets the line string to match the specified list of points.
void prepareGeometry() override
Definition: qgsgeos.cpp:158
Point geometry type, with support for z-dimension and m-values.
Definition: qgspointv2.h:34
bool isEmpty() const
QgsAbstractGeometryV2 * intersection(const QgsAbstractGeometryV2 &geom, QString *errorMsg=nullptr) const override
Definition: qgsgeos.cpp:178
const char * constData() const
static void logMessage(const QString &message, const QString &tag=QString::null, MessageLevel level=WARNING)
add a message to the instance (and create it if necessary)
virtual bool addGeometry(QgsAbstractGeometryV2 *g) override
Adds a geometry and takes ownership.
virtual QgsLineStringV2 * clone() const override
Clones the geometry by performing a deep copy.
virtual QgsPolygonV2 * toPolygon() const
QgsAbstractGeometryV2 * envelope(QString *errorMsg=nullptr) const override
Definition: qgsgeos.cpp:1358
bool disjoint(const QgsAbstractGeometryV2 &geom, QString *errorMsg=nullptr) const override
Definition: qgsgeos.cpp:277
const QgsAbstractGeometryV2 * mGeometry
QString relate(const QgsAbstractGeometryV2 &geom, QString *errorMsg=nullptr) const override
Returns the Dimensional Extended 9 Intersection Model (DE-9IM) representation of the relationship bet...
Definition: qgsgeos.cpp:282
bool touches(const QgsAbstractGeometryV2 &geom, QString *errorMsg=nullptr) const override
Definition: qgsgeos.cpp:252
QByteArray toLocal8Bit() const
void reserve(int size)
static GEOSGeometry * asGeos(const QgsAbstractGeometryV2 *geom, double precision=0)
Definition: qgsgeos.cpp:1051
virtual QString geometryType() const =0
Returns a unique string representing the geometry type.
const QgsCurveV2 * interiorRing(int i) const
void geometryChanged() override
Removes caches.
Definition: qgsgeos.cpp:149
bool crosses(const QgsAbstractGeometryV2 &geom, QString *errorMsg=nullptr) const override
Definition: qgsgeos.cpp:257
int numGeometries() const
Returns the number of geometries within the collection.
void addVertex(const QgsPointV2 &pt)
Adds a new vertex to the end of the line string.
virtual QgsLineStringV2 * curveToLine() const =0
Returns a new line string geometry corresponding to a segmentized approximation of the curve...
const_iterator constBegin() const
bool isEmpty() const
QgsAbstractGeometryV2 * reshapeGeometry(const QgsLineStringV2 &reshapeWithLine, int *errorCode, QString *errorMsg=nullptr) const
Definition: qgsgeos.cpp:1668
bool centroid(QgsPointV2 &pt, QString *errorMsg=nullptr) const override
Definition: qgsgeos.cpp:1331
QgsGeometry shortestLine(const QgsGeometry &other, QString *errorMsg=nullptr) const
Returns the shortest line joining this geometry to the other geometry.
Definition: qgsgeos.cpp:1811
virtual bool addGeometry(QgsAbstractGeometryV2 *g)
Adds a geometry and takes ownership.
virtual int dimension() const =0
Returns the inherent dimension of the geometry.
int count(const T &value) const
bool isEqual(const QgsAbstractGeometryV2 &geom, QString *errorMsg=nullptr) const override
Definition: qgsgeos.cpp:1433
const QgsAbstractGeometryV2 * geometryN(int n) const
Returns a const reference to a geometry from within the collection.
double m() const
Returns the point&#39;s m value.
Definition: qgspointv2.h:86
bool pointOnSurface(QgsPointV2 &pt, QString *errorMsg=nullptr) const override
Definition: qgsgeos.cpp:1373
double distance(const QgsAbstractGeometryV2 &geom, QString *errorMsg=nullptr) const override
Definition: qgsgeos.cpp:222
Multi polygon geometry collection.
Contains geometry relation and modification algorithms.
bool isEmpty() const
Returns true if the geometry is empty (ie, contains no underlying geometry accessible via geometry)...
Curve polygon geometry type.
bool overlaps(const QgsAbstractGeometryV2 &geom, QString *errorMsg=nullptr) const override
Definition: qgsgeos.cpp:267
QgsAbstractGeometryV2 * difference(const QgsAbstractGeometryV2 &geom, QString *errorMsg=nullptr) const override
Definition: qgsgeos.cpp:183
int numInteriorRings() const
virtual QgsAbstractGeometryV2 * clone() const =0
Clones the geometry by performing a deep copy.
Abstract base class for curved geometry type.
Definition: qgscurvev2.h:32
int size() const
QgsAbstractGeometryV2 * offsetCurve(double distance, int segments, int joinStyle, double mitreLimit, QString *errorMsg=nullptr) const override
Definition: qgsgeos.cpp:1652
#define DEFAULT_QUADRANT_SEGMENTS
Definition: qgsgeos.cpp:32
QgsAbstractGeometryV2 * convexHull(QString *errorMsg=nullptr) const override
Definition: qgsgeos.cpp:1402
QgsAbstractGeometryV2 * buffer(double distance, int segments, QString *errorMsg=nullptr) const override
Definition: qgsgeos.cpp:1263
int numPoints() const override
Returns the number of points in the curve.