[processing] port service area (from layer) alg to c++

qgis · Nov 6, 2019 · 26c90c8 · 26c90c8
1 parent 5970d1c
commit 26c90c8
Show file tree

Hide file tree

Showing 5 changed files with 334 additions and 2 deletions.
diff --git a/python/plugins/processing/algs/qgis/QgisAlgorithmProvider.py b/python/plugins/processing/algs/qgis/QgisAlgorithmProvider.py
@@ -112,7 +112,7 @@
 from .Ruggedness import Ruggedness
 from .SelectByAttribute import SelectByAttribute
 from .SelectByExpression import SelectByExpression
-from .ServiceAreaFromLayer import ServiceAreaFromLayer
+#from .ServiceAreaFromLayer import ServiceAreaFromLayer
 #from .ServiceAreaFromPoint import ServiceAreaFromPoint
 from .SetMValue import SetMValue
 from .SetRasterStyle import SetRasterStyle
@@ -221,7 +221,7 @@ def getAlgs(self):
                 Ruggedness(),
                 SelectByAttribute(),
                 SelectByExpression(),
-                ServiceAreaFromLayer(),
+                #ServiceAreaFromLayer(),
                 #ServiceAreaFromPoint(),
                 SetMValue(),
                 SetRasterStyle(),

diff --git a/src/analysis/CMakeLists.txt b/src/analysis/CMakeLists.txt
@@ -100,6 +100,7 @@ SET(QGIS_ANALYSIS_SRCS
   processing/qgsalgorithmrotate.cpp
   processing/qgsalgorithmsaveselectedfeatures.cpp
   processing/qgsalgorithmsegmentize.cpp
+  processing/qgsalgorithmserviceareafromlayer.cpp
   processing/qgsalgorithmserviceareafrompoint.cpp
   processing/qgsalgorithmshortestpathlayertopoint.cpp
   processing/qgsalgorithmshortestpathpointtolayer.cpp

diff --git a/src/analysis/processing/qgsalgorithmserviceareafromlayer.cpp b/src/analysis/processing/qgsalgorithmserviceareafromlayer.cpp
@@ -0,0 +1,276 @@
+/***************************************************************************
+                         qgsalgorithmserviceareafromlayer.cpp
+                         ---------------------
+    begin                : July 2018
+    copyright            : (C) 2018 by Alexander Bruy
+    email                : alexander dot bruy 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 "qgsalgorithmserviceareafromlayer.h"
+
+#include "qgsgeometryutils.h"
+#include "qgsgraphanalyzer.h"
+
+///@cond PRIVATE
+
+QString QgsServiceAreaFromLayerAlgorithm::name() const
+{
+  return QStringLiteral( "serviceareafromlayer" );
+}
+
+QString QgsServiceAreaFromLayerAlgorithm::displayName() const
+{
+  return QObject::tr( "Service area (from layer)" );
+}
+
+QStringList QgsServiceAreaFromLayerAlgorithm::tags() const
+{
+  return QObject::tr( "network,service,area,shortest,fastest" ).split( ',' );
+}
+
+QString QgsServiceAreaFromLayerAlgorithm::shortHelpString() const
+{
+  return QObject::tr( "This algorithm computes service areas from the set of start points defined by point layer using travel time or distance." );
+}
+
+QgsServiceAreaFromLayerAlgorithm *QgsServiceAreaFromLayerAlgorithm::createInstance() const
+{
+  return new QgsServiceAreaFromLayerAlgorithm();
+}
+
+void QgsServiceAreaFromLayerAlgorithm::initAlgorithm( const QVariantMap & )
+{
+  addCommonParams();
+  addParameter( new QgsProcessingParameterFeatureSource( QStringLiteral( "START_POINTS" ), QObject::tr( "Vector layer with start points" ), QList< int >() << QgsProcessing::TypeVectorPoint ) );
+  addParameter( new QgsProcessingParameterNumber( QStringLiteral( "TRAVEL_COST" ), QObject::tr( "Travel cost (distance for 'Shortest', time for 'Fastest')" ),
+                QgsProcessingParameterNumber::Double, 0, false, 0 ) );
+
+  std::unique_ptr< QgsProcessingParameterBoolean > includeBounds = qgis::make_unique< QgsProcessingParameterBoolean >( QStringLiteral( "INCLUDE_BOUNDS" ), QObject::tr( "Include upper/lower bound points" ), false, true );
+  includeBounds->setFlags( includeBounds->flags() | QgsProcessingParameterDefinition::FlagAdvanced );
+  addParameter( includeBounds.release() );
+
+  std::unique_ptr< QgsProcessingParameterFeatureSink > outputLines = qgis::make_unique< QgsProcessingParameterFeatureSink >( QStringLiteral( "OUTPUT_LINES" ),  QObject::tr( "Service area (lines)" ),
+      QgsProcessing::TypeVectorLine, QVariant(), true );
+  outputLines->setCreateByDefault( true );
+  addParameter( outputLines.release() );
+
+  std::unique_ptr< QgsProcessingParameterFeatureSink > outputPoints = qgis::make_unique< QgsProcessingParameterFeatureSink >( QStringLiteral( "OUTPUT" ),  QObject::tr( "Service area (boundary nodes)" ),
+      QgsProcessing::TypeVectorPoint, QVariant(), true );
+  outputPoints->setCreateByDefault( false );
+  addParameter( outputPoints.release() );
+}
+
+QVariantMap QgsServiceAreaFromLayerAlgorithm::processAlgorithm( const QVariantMap &parameters, QgsProcessingContext &context, QgsProcessingFeedback *feedback )
+{
+  loadCommonParams( parameters, context, feedback );
+
+  std::unique_ptr< QgsFeatureSource > startPoints( parameterAsSource( parameters, QStringLiteral( "START_POINTS" ), context ) );
+  if ( !startPoints )
+    throw QgsProcessingException( invalidSourceError( parameters, QStringLiteral( "START_POINTS" ) ) );
+
+  double travelCost = parameterAsDouble( parameters, QStringLiteral( "TRAVEL_COST" ), context );
+
+  bool includeBounds = true;  // default to true to maintain 3.0 API
+  if ( parameters.contains( QStringLiteral( "INCLUDE_BOUNDS" ) ) )
+  {
+    includeBounds = parameterAsBool( parameters, QStringLiteral( "INCLUDE_BOUNDS" ), context );
+  }
+
+  QVector< QgsPointXY > points;
+  QHash< int, QgsAttributes > sourceAttributes;
+  loadPoints( startPoints.get(), points, sourceAttributes, context, feedback );
+
+  feedback->pushInfo( QObject::tr( "Building graph…" ) );
+  QVector< QgsPointXY > snappedPoints;
+  mDirector->makeGraph( mBuilder.get(), points, snappedPoints, feedback );
+
+  feedback->pushInfo( QObject::tr( "Calculating service areas…" ) );
+  QgsGraph *graph = mBuilder->graph();
+
+  QgsFields fields = startPoints->fields();
+  fields.append( QgsField( QStringLiteral( "type" ), QVariant::String ) );
+  fields.append( QgsField( QStringLiteral( "start" ), QVariant::String ) );
+
+  QString pointsSinkId;
+  std::unique_ptr< QgsFeatureSink > pointsSink( parameterAsSink( parameters, QStringLiteral( "OUTPUT" ), context, pointsSinkId, fields,
+      QgsWkbTypes::MultiPoint, mNetwork->sourceCrs() ) );
+
+  QString linesSinkId;
+  std::unique_ptr< QgsFeatureSink > linesSink( parameterAsSink( parameters, QStringLiteral( "OUTPUT_LINES" ), context, linesSinkId, fields,
+      QgsWkbTypes::MultiLineString, mNetwork->sourceCrs() ) );
+
+  int idxStart;
+  QString origPoint;
+  QVector< int > tree;
+  QVector< double > costs;
+
+  int inboundEdgeIndex;
+  double startVertexCost, endVertexCost;
+  QgsPointXY startPoint, endPoint;
+  QgsGraphEdge edge;
+
+  QgsFeature feat;
+  QgsAttributes attributes;
+
+  int step =  snappedPoints.size() > 0 ? 100.0 / snappedPoints.size() : 1;
+  for ( int i = 0; i < snappedPoints.size(); i++ )
+  {
+    if ( feedback->isCanceled() )
+    {
+      break;
+    }
+
+    idxStart = graph->findVertex( snappedPoints.at( i ) );
+    origPoint = points.at( i ).toString();
+
+    QgsGraphAnalyzer::dijkstra( graph, idxStart, 0, &tree, &costs );
+
+    QgsMultiPointXY areaPoints;
+    QgsMultiPolylineXY lines;
+    QSet< int > vertices;
+
+    for ( int j = 0; j < costs.size(); j++ )
+    {
+      inboundEdgeIndex = tree.at( j );
+
+      if ( inboundEdgeIndex == -1 and j != idxStart )
+      {
+        // unreachable vertex
+        continue;
+      }
+
+      startVertexCost = costs.at( j );
+      if ( startVertexCost > travelCost )
+      {
+        // vertex is too expensive, discard
+        continue;
+      }
+
+      vertices.insert( j );
+      startPoint = graph->vertex( j ).point();
+
+      // find all edges coming from this vertex
+      for ( int edgeId : graph->vertex( j ).outgoingEdges() )
+      {
+        edge = graph->edge( edgeId );
+        endVertexCost = startVertexCost + edge.cost( 0 ).toDouble();
+        endPoint = graph->vertex( edge.toVertex() ).point();
+        if ( endVertexCost <= travelCost )
+        {
+          // end vertex is cheap enough to include
+          vertices.insert( edge.toVertex() );
+          lines.push_back( QgsPolylineXY() << startPoint << endPoint );
+        }
+        else
+        {
+          // travelCost sits somewhere on this edge, interpolate position
+          QgsPointXY interpolatedEndPoint = QgsGeometryUtils::interpolatePointOnLineByValue( startPoint.x(), startPoint.y(), startVertexCost,
+                                            endPoint.x(), endPoint.y(), endVertexCost, travelCost );
+
+          areaPoints.push_back( interpolatedEndPoint );
+          lines.push_back( QgsPolylineXY() << startPoint << interpolatedEndPoint );
+        }
+      } // edges
+    } // costs
+
+    // convert to list so he have same order of points between algorithm runs
+    QList< int > verticesList = vertices.toList();
+    std::sort( verticesList.begin(), verticesList.end() );
+    for ( int v : verticesList )
+    {
+      areaPoints.push_back( graph->vertex( v ).point() );
+    }
+    // this produces correct results, but order of point is not the same, so tests failed
+    //QSet< int >::const_iterator it = vertices.constBegin();
+    //while ( it != vertices.constEnd() )
+    //{
+    //  areaPoints.push_back( graph->vertex( *it ).point() );
+    //  it++;
+    //}
+
+    if ( pointsSink )
+    {
+      QgsGeometry geomPoints = QgsGeometry::fromMultiPointXY( areaPoints );
+      feat.setGeometry( geomPoints );
+      attributes = sourceAttributes.value( i + 1 );
+      attributes << QStringLiteral( "within" ) << origPoint;
+      feat.setAttributes( attributes );
+      pointsSink->addFeature( feat, QgsFeatureSink::FastInsert );
+
+      if ( includeBounds )
+      {
+        QgsMultiPointXY upperBoundary, lowerBoundary;
+        QVector< int > nodes;
+
+        int vertexId;
+        for ( int v = 0; v < costs.size(); v++ )
+        {
+          if ( costs.at( v ) > travelCost and tree.at( v ) != -1 )
+          {
+            vertexId = graph->edge( tree.at( v ) ).fromVertex();
+            if ( costs.at( vertexId ) <= travelCost )
+            {
+              nodes.push_back( v );
+            }
+          }
+        } // costs
+
+        for ( int n : nodes )
+        {
+          upperBoundary.push_back( graph->vertex( graph->edge( tree.at( n ) ).toVertex() ).point() );
+          lowerBoundary.push_back( graph->vertex( graph->edge( tree.at( n ) ).fromVertex() ).point() );
+        } // nodes
+
+        QgsGeometry geomUpper = QgsGeometry::fromMultiPointXY( upperBoundary );
+        QgsGeometry geomLower = QgsGeometry::fromMultiPointXY( lowerBoundary );
+
+        feat.setGeometry( geomUpper );
+        attributes = sourceAttributes.value( i + 1 );
+        attributes << QStringLiteral( "upper" ) << origPoint;
+        feat.setAttributes( attributes );
+        pointsSink->addFeature( feat, QgsFeatureSink::FastInsert );
+
+        feat.setGeometry( geomLower );
+        attributes = sourceAttributes.value( i + 1 );
+        attributes << QStringLiteral( "lower" ) << origPoint;
+        feat.setAttributes( attributes );
+        pointsSink->addFeature( feat, QgsFeatureSink::FastInsert );
+      } // includeBounds
+    }
+
+    if ( linesSink )
+    {
+      QgsGeometry geomLines = QgsGeometry::fromMultiPolylineXY( lines );
+      feat.setGeometry( geomLines );
+      attributes = sourceAttributes.value( i + 1 );
+      attributes << QStringLiteral( "lines" ) << origPoint;
+      feat.setAttributes( attributes );
+      linesSink->addFeature( feat, QgsFeatureSink::FastInsert );
+    }
+
+    feedback->setProgress( i * step );
+  } // snappedPoints
+
+  QVariantMap outputs;
+  if ( pointsSink )
+  {
+    outputs.insert( QStringLiteral( "OUTPUT" ), pointsSinkId );
+  }
+  if ( linesSink )
+  {
+    outputs.insert( QStringLiteral( "OUTPUT_LINES" ), linesSinkId );
+  }
+
+  return outputs;
+}
+
+///@endcond
diff --git a/src/analysis/processing/qgsalgorithmserviceareafromlayer.h b/src/analysis/processing/qgsalgorithmserviceareafromlayer.h
@@ -0,0 +1,53 @@
+/***************************************************************************
+                         qgsalgorithmserviceareafromlayer.h
+                         ---------------------
+    begin                : JUly 2018
+    copyright            : (C) 2018 by Alexander Bruy
+    email                : alexander dot bruy 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 QGSALGORITHMSERVICEAREAFROMLAYER_H
+#define QGSALGORITHMSERVICEAREAFROMLAYER_H
+
+#define SIP_NO_FILE
+
+#include "qgis.h"
+#include "qgsalgorithmnetworkanalysisbase.h"
+
+///@cond PRIVATE
+
+/**
+ * Native service area (from layer) algorithm.
+ */
+class QgsServiceAreaFromLayerAlgorithm : public QgsNetworkAnalysisAlgorithmBase
+{
+
+  public:
+
+    QgsServiceAreaFromLayerAlgorithm() = default;
+    void initAlgorithm( const QVariantMap &configuration = QVariantMap() ) override;
+    QString name() const override;
+    QString displayName() const override;
+    QStringList tags() const override;
+    QString shortHelpString() const override;
+    QgsServiceAreaFromLayerAlgorithm *createInstance() const override SIP_FACTORY;
+
+  protected:
+
+    QVariantMap processAlgorithm( const QVariantMap &parameters,
+                                  QgsProcessingContext &context, QgsProcessingFeedback *feedback ) override;
+
+};
+
+///@endcond PRIVATE
+
+#endif // QGSALGORITHMSERVICEAREAFROMLAYER_H
diff --git a/src/analysis/processing/qgsnativealgorithms.cpp b/src/analysis/processing/qgsnativealgorithms.cpp
@@ -94,6 +94,7 @@
 #include "qgsalgorithmrotate.h"
 #include "qgsalgorithmsaveselectedfeatures.h"
 #include "qgsalgorithmsegmentize.h"
+#include "qgsalgorithmserviceareafromlayer.h"
 #include "qgsalgorithmserviceareafrompoint.h"
 #include "qgsalgorithmshortestpathlayertopoint.h"
 #include "qgsalgorithmshortestpathpointtolayer.h"
@@ -251,6 +252,7 @@ void QgsNativeAlgorithms::loadAlgorithms()
   addAlgorithm( new QgsSegmentizeByMaximumAngleAlgorithm() );
   addAlgorithm( new QgsSegmentizeByMaximumDistanceAlgorithm() );
   addAlgorithm( new QgsSelectByLocationAlgorithm() );
+  addAlgorithm( new QgsServiceAreaFromLayerAlgorithm() );
   addAlgorithm( new QgsServiceAreaFromPointAlgorithm() );
   addAlgorithm( new QgsShortestPathLayerToPointAlgorithm() );
   addAlgorithm( new QgsShortestPathPointToLayerAlgorithm() );