Skip to content

Commit

Permalink
Expose GEOS coverage simplify as a processing algorithm
Browse files Browse the repository at this point in the history
  • Loading branch information
nyalldawson committed Oct 15, 2023
1 parent e158f43 commit 6671d42
Show file tree
Hide file tree
Showing 7 changed files with 353 additions and 0 deletions.
46 changes: 46 additions & 0 deletions python/core/auto_generated/geometry/qgsgeometry.sip.in
Original file line number Diff line number Diff line change
Expand Up @@ -1808,6 +1808,52 @@ instead of polygons.
An empty geometry will be returned if the diagram could not be calculated.

.. versionadded:: 3.0
%End

Qgis::CoverageValidityResult validateCoverage( double gapWidth, QgsGeometry *invalidEdges /Out/ = 0 ) const throw( QgsNotSupportedException );
%Docstring
Analyze a coverage (represented as a collection of polygonal geometry with exactly matching edge
geometry) to find places where the assumption of exactly matching edges is not met.

The input geometry is the polygonal coverage to access, stored in a geometry collection. All members must be POLYGON
or MULTIPOLYGON.

:param gapWidth: The maximum width of gaps to detect.

:return: - validity check result
- invalidEdges: When there are invalidities in the coverage, will be set with a geometry collection of the same length as the input, with a MULTILINESTRING of the error edges for each invalid polygon, or an EMPTY where the polygon is a valid participant in the coverage.

This method requires a QGIS build based on GEOS 3.12 or later.

:raises QgsNotSupportedException: on QGIS builds based on GEOS 3.11 or earlier.

.. seealso:: :py:func:`simplifyCoverageVW`

.. versionadded:: 3.36
%End

QgsGeometry simplifyCoverageVW( double tolerance, bool preserveBoundary ) const throw( QgsNotSupportedException );
%Docstring
Operates on a coverage (represented as a list of polygonal geometry with exactly matching edge geometry) to apply
a Visvalingam–Whyatt simplification to the edges, reducing complexity in proportion with the provided tolerance,
while retaining a valid coverage (no edges will cross or touch after the simplification).

Geometries never disappear, but they may be simplified down to just a triangle. Also, some invalid geoms
(such as Polygons which have too few non-repeated points) will be returned unchanged.

If the input dataset is not a valid coverage due to overlaps, it will still be simplified, but invalid topology
such as crossing edges will still be invalid.

:param tolerance: A tolerance parameter in linear units.
:param preserveBoundary: Set to ``True`` to preserve the outside edges of the coverage without simplification, or ``False`` to allow them to be simplified.

This method requires a QGIS build based on GEOS 3.12 or later.

:raises QgsNotSupportedException: on QGIS builds based on GEOS 3.11 or earlier.

.. seealso:: :py:func:`validateCoverage`

.. versionadded:: 3.36
%End

QgsGeometry node() const;
Expand Down
1 change: 1 addition & 0 deletions src/analysis/CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -66,6 +66,7 @@ set(QGIS_ANALYSIS_SRCS
processing/qgsalgorithmconstantraster.cpp
processing/qgsalgorithmconverttocurves.cpp
processing/qgsalgorithmconvexhull.cpp
processing/qgsalgorithmcoveragesimplify.cpp
processing/qgsalgorithmcreatedirectory.cpp
processing/qgsalgorithmdbscanclustering.cpp
processing/qgsalgorithmdelaunaytriangulation.cpp
Expand Down
175 changes: 175 additions & 0 deletions src/analysis/processing/qgsalgorithmcoveragesimplify.cpp
Original file line number Diff line number Diff line change
@@ -0,0 +1,175 @@
/***************************************************************************
qgsalgorithmcoveragesimplify.cpp
---------------------
begin : October 2023
copyright : (C) 2023 by Nyall Dawson
email : nyall dot dawson at gmail dot com
***************************************************************************/

/***************************************************************************
* *
* This program is free software; you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by *
* the Free Software Foundation; either version 2 of the License, or *
* (at your option) any later version. *
* *
***************************************************************************/


#include "qgsalgorithmcoveragesimplify.h"
#include "qgsgeometrycollection.h"
#include "qgsgeos.h"


///@cond PRIVATE

QString QgsCoverageSimplifyAlgorithm::name() const
{
return QStringLiteral( "coveragesimplify" );
}

QString QgsCoverageSimplifyAlgorithm::displayName() const
{
return QObject::tr( "Simplify coverage" );
}

QStringList QgsCoverageSimplifyAlgorithm::tags() const
{
return QObject::tr( "topological,boundary" ).split( ',' );
}

QString QgsCoverageSimplifyAlgorithm::group() const
{
return QObject::tr( "Vector geometry" );
}

QString QgsCoverageSimplifyAlgorithm::groupId() const
{
return QStringLiteral( "vectorgeometry" );
}

void QgsCoverageSimplifyAlgorithm::initAlgorithm( const QVariantMap & )
{
addParameter( new QgsProcessingParameterFeatureSource( QStringLiteral( "INPUT" ), QObject::tr( "Input layer" ), QList< int >() << QgsProcessing::TypeVectorPolygon ) );
addParameter( new QgsProcessingParameterDistance( QStringLiteral( "TOLERANCE" ),
QObject::tr( "Tolerance" ), 1.0, QStringLiteral( "INPUT" ), false, 0, 10000000.0) );
addParameter( new QgsProcessingParameterBoolean( QStringLiteral( "PRESERVE_BOUNDARY" ), QObject::tr( "Preserve boundary" ), false ) );

addParameter( new QgsProcessingParameterFeatureSink( QStringLiteral( "OUTPUT" ), QObject::tr( "Simplified" ) ) );
}

QString QgsCoverageSimplifyAlgorithm::shortHelpString() const
{
return QObject::tr( "" );
}

QgsCoverageSimplifyAlgorithm *QgsCoverageSimplifyAlgorithm::createInstance() const
{
return new QgsCoverageSimplifyAlgorithm();
}

QVariantMap QgsCoverageSimplifyAlgorithm::processAlgorithm( const QVariantMap &parameters, QgsProcessingContext &context, QgsProcessingFeedback *feedback )
{
std::unique_ptr< QgsProcessingFeatureSource > source( parameterAsSource( parameters, QStringLiteral( "INPUT" ), context ) );
if ( !source )
throw QgsProcessingException( invalidSourceError( parameters, QStringLiteral( "INPUT" ) ) );

const bool preserveBoundary = parameterAsBoolean( parameters, QStringLiteral( "PRESERVE_BOUNDARY" ), context );
const double tolerance = parameterAsDouble( parameters, QStringLiteral( "TOLERANCE" ), context );

QString sinkId;
std::unique_ptr< QgsFeatureSink > sink( parameterAsSink( parameters, QStringLiteral( "OUTPUT" ), context, sinkId, source->fields(), source->wkbType(), source->sourceCrs() ) );
if ( !sink )
throw QgsProcessingException( invalidSinkError( parameters, QStringLiteral( "OUTPUT" ) ) );

// we have no choice but to build up a list of features in advance
QVector< QgsFeature > featuresWithGeom;
QVector< QgsFeature > featuresWithoutGeom;
QgsGeometryCollection collection;

const long count = source->featureCount();
if ( count > 0 )
{
featuresWithGeom.reserve( count );
collection.reserve( count );
}

const double step = count > 0 ? 100.0 / count : 1;
int current = 0;

feedback->pushInfo( QObject::tr("Collecting features"));

QgsFeature inFeature;
QgsFeatureIterator features = source->getFeatures();
while ( features.nextFeature( inFeature ) )
{
if ( feedback->isCanceled() )
{
break;
}

if ( inFeature.hasGeometry() )
{
featuresWithGeom.append( inFeature );
collection.addGeometry( inFeature.geometry().constGet()->clone() );
}
else
{
featuresWithoutGeom.append( inFeature );
}


feedback->setProgress( current * step * 0.2 );
current++;
}

feedback->pushInfo( QObject::tr("Simplifying coverage"));

QgsGeos geos( &collection );
QString error;
std::unique_ptr< QgsAbstractGeometry > simplified;
try
{
simplified.reset( geos.simplifyCoverageVW( tolerance, preserveBoundary, &error ) );
}
catch ( QgsNotSupportedException& e )
{
throw QgsProcessingException( e.what() );
}

if ( !simplified )
{
if ( !error.isEmpty())
throw QgsProcessingException( error );
else
throw QgsProcessingException( QObject::tr("No geometry was returned for simplified coverage") );
}

feedback->setProgress( 80 );

feedback->pushInfo( QObject::tr("Storing features"));
long long featureIndex = 0;
for ( auto partsIt = simplified->const_parts_begin(); partsIt != simplified->const_parts_end(); ++partsIt )
{
QgsFeature outFeature = featuresWithGeom.value( featureIndex );
outFeature.setGeometry( QgsGeometry( *partsIt ? (*partsIt )->clone() : nullptr ));
if ( !sink->addFeature( outFeature, QgsFeatureSink::FastInsert ) )
throw QgsProcessingException( writeFeatureError( sink.get(), parameters, QStringLiteral( "OUTPUT" ) ) );


feedback->setProgress( featureIndex * step * 0.2 + 80 );
featureIndex++;
}
for ( const QgsFeature& feature : featuresWithoutGeom)
{
QgsFeature outFeature =feature;
if ( !sink->addFeature( outFeature, QgsFeatureSink::FastInsert ) )
throw QgsProcessingException( writeFeatureError( sink.get(), parameters, QStringLiteral( "OUTPUT" ) ) );
}

QVariantMap outputs;
outputs.insert( QStringLiteral( "OUTPUT" ), sinkId );
return outputs;
}

///@endcond
55 changes: 55 additions & 0 deletions src/analysis/processing/qgsalgorithmcoveragesimplify.h
Original file line number Diff line number Diff line change
@@ -0,0 +1,55 @@
/***************************************************************************
qgsalgorithmcoveragesimplify.h
---------------------
begin : October 2023
copyright : (C) 2023 by Nyall Dawson
email : nyall dot dawson at gmail dot com
***************************************************************************/

/***************************************************************************
* *
* This program is free software; you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by *
* the Free Software Foundation; either version 2 of the License, or *
* (at your option) any later version. *
* *
***************************************************************************/

#ifndef QGSALGORITHMCOVERAGESIMPLIFY_H
#define QGSALGORITHMCOVERAGESIMPLIFY_H

#define SIP_NO_FILE

#include "qgis_sip.h"
#include "qgsprocessingalgorithm.h"

///@cond PRIVATE

/**
* Native coverage simplification algorithm.
*/
class QgsCoverageSimplifyAlgorithm : public QgsProcessingAlgorithm
{
public:
QgsCoverageSimplifyAlgorithm() = default;
void initAlgorithm( const QVariantMap &configuration = QVariantMap() ) override;
QString name() const override;
QString displayName() const override;
QStringList tags() const override;
QString group() const override;
QString groupId() const override;
QString shortHelpString() const override;
QgsCoverageSimplifyAlgorithm *createInstance() const override SIP_FACTORY;

protected:

QVariantMap processAlgorithm( const QVariantMap &parameters,
QgsProcessingContext &context, QgsProcessingFeedback *feedback ) override;
};


///@endcond PRIVATE

#endif // QGSALGORITHMCOVERAGESIMPLIFY_H


2 changes: 2 additions & 0 deletions src/analysis/processing/qgsnativealgorithms.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -47,6 +47,7 @@
#include "qgsalgorithmconstantraster.h"
#include "qgsalgorithmconverttocurves.h"
#include "qgsalgorithmconvexhull.h"
#include "qgsalgorithmcoveragesimplify.h"
#include "qgsalgorithmcreatedirectory.h"
#include "qgsalgorithmdbscanclustering.h"
#include "qgsalgorithmdelaunaytriangulation.h"
Expand Down Expand Up @@ -313,6 +314,7 @@ void QgsNativeAlgorithms::loadAlgorithms()
addAlgorithm( new QgsConstantRasterAlgorithm() );
addAlgorithm( new QgsConvertToCurvesAlgorithm() );
addAlgorithm( new QgsConvexHullAlgorithm() );
addAlgorithm( new QgsCoverageSimplifyAlgorithm() );
addAlgorithm( new QgsCreateDirectoryAlgorithm() );
addAlgorithm( new QgsDbscanClusteringAlgorithm() );
addAlgorithm( new QgsDelaunayTriangulationAlgorithm() );
Expand Down
33 changes: 33 additions & 0 deletions src/core/geometry/qgsgeometry.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -2638,6 +2638,39 @@ QgsGeometry QgsGeometry::delaunayTriangulation( double tolerance, bool edgesOnly
return result;
}

Qgis::CoverageValidityResult QgsGeometry::validateCoverage( double gapWidth, QgsGeometry *invalidEdges ) const
{
if ( !d->geometry )
{
return Qgis::CoverageValidityResult::Error;
}

QgsGeos geos( d->geometry.get() );
mLastError.clear();
std::unique_ptr< QgsAbstractGeometry > invalidEdgesGeom;

const Qgis::CoverageValidityResult result = geos.validateCoverage( gapWidth, invalidEdges ? &invalidEdgesGeom : nullptr, &mLastError );

if ( invalidEdges && invalidEdgesGeom )
*invalidEdges = QgsGeometry( std::move( invalidEdgesGeom ) );

return result;
}

QgsGeometry QgsGeometry::simplifyCoverageVW( double tolerance, bool preserveBoundary ) const
{
if ( !d->geometry )
{
return QgsGeometry();
}

QgsGeos geos( d->geometry.get() );
mLastError.clear();
QgsGeometry result( geos.simplifyCoverageVW( tolerance, preserveBoundary, &mLastError ) );
result.mLastError = mLastError;
return result;
}

QgsGeometry QgsGeometry::node() const
{
if ( !d->geometry )
Expand Down
41 changes: 41 additions & 0 deletions src/core/geometry/qgsgeometry.h
Original file line number Diff line number Diff line change
Expand Up @@ -1833,6 +1833,47 @@ class CORE_EXPORT QgsGeometry
*/
QgsGeometry delaunayTriangulation( double tolerance = 0.0, bool edgesOnly = false ) const;

/**
* Analyze a coverage (represented as a collection of polygonal geometry with exactly matching edge
* geometry) to find places where the assumption of exactly matching edges is not met.
*
* The input geometry is the polygonal coverage to access, stored in a geometry collection. All members must be POLYGON
* or MULTIPOLYGON.
*
* \param gapWidth The maximum width of gaps to detect.
* \param invalidEdges When there are invalidities in the coverage, will be set with a geometry collection of the same length as the input, with a MULTILINESTRING of the error edges for each invalid polygon, or an EMPTY where the polygon is a valid participant in the coverage.
* \returns validity check result
*
* This method requires a QGIS build based on GEOS 3.12 or later.
*
* \throws QgsNotSupportedException on QGIS builds based on GEOS 3.11 or earlier.
* \see simplifyCoverageVW()
* \since QGIS 3.36
*/
Qgis::CoverageValidityResult validateCoverage( double gapWidth, QgsGeometry *invalidEdges SIP_OUT = nullptr ) const SIP_THROW( QgsNotSupportedException );

/**
* Operates on a coverage (represented as a list of polygonal geometry with exactly matching edge geometry) to apply
* a Visvalingam–Whyatt simplification to the edges, reducing complexity in proportion with the provided tolerance,
* while retaining a valid coverage (no edges will cross or touch after the simplification).
*
* Geometries never disappear, but they may be simplified down to just a triangle. Also, some invalid geoms
* (such as Polygons which have too few non-repeated points) will be returned unchanged.
*
* If the input dataset is not a valid coverage due to overlaps, it will still be simplified, but invalid topology
* such as crossing edges will still be invalid.
*
* \param tolerance A tolerance parameter in linear units.
* \param preserveBoundary Set to TRUE to preserve the outside edges of the coverage without simplification, or FALSE to allow them to be simplified.
*
* This method requires a QGIS build based on GEOS 3.12 or later.
*
* \throws QgsNotSupportedException on QGIS builds based on GEOS 3.11 or earlier.
* \see validateCoverage()
* \since QGIS 3.36
*/
QgsGeometry simplifyCoverageVW( double tolerance, bool preserveBoundary ) const SIP_THROW( QgsNotSupportedException );

/**
* Returns a (Multi)LineString representing the fully noded version of a collection of linestrings.
*
Expand Down

0 comments on commit 6671d42

Please sign in to comment.