QGIS API Documentation 3.28.0-Firenze (ed3ad0430f)
qgsalgorithmserviceareafromlayer.cpp
Go to the documentation of this file.
1/***************************************************************************
2 qgsalgorithmserviceareafromlayer.cpp
3 ---------------------
4 begin : July 2018
5 copyright : (C) 2018 by Alexander Bruy
6 email : alexander dot bruy at gmail dot com
7 ***************************************************************************/
8
9/***************************************************************************
10 * *
11 * This program is free software; you can redistribute it and/or modify *
12 * it under the terms of the GNU General Public License as published by *
13 * the Free Software Foundation; either version 2 of the License, or *
14 * (at your option) any later version. *
15 * *
16 ***************************************************************************/
17
19
20#include "qgsgeometryutils.h"
21#include "qgsgraphanalyzer.h"
22
24
25QString QgsServiceAreaFromLayerAlgorithm::name() const
26{
27 return QStringLiteral( "serviceareafromlayer" );
28}
29
30QString QgsServiceAreaFromLayerAlgorithm::displayName() const
31{
32 return QObject::tr( "Service area (from layer)" );
33}
34
35QStringList QgsServiceAreaFromLayerAlgorithm::tags() const
36{
37 return QObject::tr( "network,service,area,shortest,fastest" ).split( ',' );
38}
39
40QString QgsServiceAreaFromLayerAlgorithm::shortHelpString() const
41{
42 return QObject::tr( "This algorithm creates a new vector with all the edges or parts of "
43 "edges of a network line layer that can be reached within a distance "
44 "or a time, starting from features of a point layer. The distance and "
45 "the time (both referred to as \"travel cost\") must be specified "
46 "respectively in the network layer units or in hours." );
47}
48
49QgsServiceAreaFromLayerAlgorithm *QgsServiceAreaFromLayerAlgorithm::createInstance() const
50{
51 return new QgsServiceAreaFromLayerAlgorithm();
52}
53
54void QgsServiceAreaFromLayerAlgorithm::initAlgorithm( const QVariantMap & )
55{
56 addCommonParams();
57 addParameter( new QgsProcessingParameterFeatureSource( QStringLiteral( "START_POINTS" ), QObject::tr( "Vector layer with start points" ), QList< int >() << QgsProcessing::TypeVectorPoint ) );
58
59 std::unique_ptr< QgsProcessingParameterNumber > travelCost = std::make_unique< QgsProcessingParameterNumber >( QStringLiteral( "TRAVEL_COST" ), QObject::tr( "Travel cost (distance for 'Shortest', time for 'Fastest')" ), QgsProcessingParameterNumber::Double, 0, true, 0 );
60 travelCost->setFlags( travelCost->flags() | QgsProcessingParameterDefinition::FlagHidden );
61 addParameter( travelCost.release() );
62
63 addParameter( new QgsProcessingParameterNumber( QStringLiteral( "TRAVEL_COST2" ), QObject::tr( "Travel cost (distance for 'Shortest', time for 'Fastest')" ),
65
66 std::unique_ptr< QgsProcessingParameterBoolean > includeBounds = std::make_unique< QgsProcessingParameterBoolean >( QStringLiteral( "INCLUDE_BOUNDS" ), QObject::tr( "Include upper/lower bound points" ), false, true );
67 includeBounds->setFlags( includeBounds->flags() | QgsProcessingParameterDefinition::FlagAdvanced );
68 addParameter( includeBounds.release() );
69
70 std::unique_ptr< QgsProcessingParameterFeatureSink > outputLines = std::make_unique< QgsProcessingParameterFeatureSink >( QStringLiteral( "OUTPUT_LINES" ), QObject::tr( "Service area (lines)" ),
71 QgsProcessing::TypeVectorLine, QVariant(), true );
72 outputLines->setCreateByDefault( true );
73 addParameter( outputLines.release() );
74
75 std::unique_ptr< QgsProcessingParameterFeatureSink > outputPoints = std::make_unique< QgsProcessingParameterFeatureSink >( QStringLiteral( "OUTPUT" ), QObject::tr( "Service area (boundary nodes)" ),
76 QgsProcessing::TypeVectorPoint, QVariant(), true );
77 outputPoints->setCreateByDefault( false );
78 addParameter( outputPoints.release() );
79}
80
81QVariantMap QgsServiceAreaFromLayerAlgorithm::processAlgorithm( const QVariantMap &parameters, QgsProcessingContext &context, QgsProcessingFeedback *feedback )
82{
83 loadCommonParams( parameters, context, feedback );
84
85 std::unique_ptr< QgsFeatureSource > startPoints( parameterAsSource( parameters, QStringLiteral( "START_POINTS" ), context ) );
86 if ( !startPoints )
87 throw QgsProcessingException( invalidSourceError( parameters, QStringLiteral( "START_POINTS" ) ) );
88
89 // use older deprecated travel cost style if specified, to maintain old api
90 const bool useOldTravelCost = parameters.value( QStringLiteral( "TRAVEL_COST" ) ).isValid();
91 double travelCost = parameterAsDouble( parameters, useOldTravelCost ? QStringLiteral( "TRAVEL_COST" ) : QStringLiteral( "TRAVEL_COST2" ), context );
92
93 int strategy = parameterAsInt( parameters, QStringLiteral( "STRATEGY" ), context );
94 if ( strategy && !useOldTravelCost )
95 travelCost *= mMultiplier;
96
97 bool includeBounds = true; // default to true to maintain 3.0 API
98 if ( parameters.contains( QStringLiteral( "INCLUDE_BOUNDS" ) ) )
99 {
100 includeBounds = parameterAsBool( parameters, QStringLiteral( "INCLUDE_BOUNDS" ), context );
101 }
102
103 QVector< QgsPointXY > points;
104 QHash< int, QgsAttributes > sourceAttributes;
105 loadPoints( startPoints.get(), points, sourceAttributes, context, feedback );
106
107 feedback->pushInfo( QObject::tr( "Building graph…" ) );
108 QVector< QgsPointXY > snappedPoints;
109 mDirector->makeGraph( mBuilder.get(), points, snappedPoints, feedback );
110
111 feedback->pushInfo( QObject::tr( "Calculating service areas…" ) );
112 std::unique_ptr< QgsGraph > graph( mBuilder->takeGraph() );
113
114 QgsFields fields = startPoints->fields();
115 fields.append( QgsField( QStringLiteral( "type" ), QVariant::String ) );
116 fields.append( QgsField( QStringLiteral( "start" ), QVariant::String ) );
117
118 QString pointsSinkId;
119 std::unique_ptr< QgsFeatureSink > pointsSink( parameterAsSink( parameters, QStringLiteral( "OUTPUT" ), context, pointsSinkId, fields,
120 QgsWkbTypes::MultiPoint, mNetwork->sourceCrs() ) );
121
122 QString linesSinkId;
123 std::unique_ptr< QgsFeatureSink > linesSink( parameterAsSink( parameters, QStringLiteral( "OUTPUT_LINES" ), context, linesSinkId, fields,
124 QgsWkbTypes::MultiLineString, mNetwork->sourceCrs() ) );
125
126 int idxStart;
127 QString origPoint;
128 QVector< int > tree;
129 QVector< double > costs;
130
131 int inboundEdgeIndex;
132 double startVertexCost, endVertexCost;
133 QgsPointXY startPoint, endPoint;
134 QgsGraphEdge edge;
135
136 QgsFeature feat;
137 QgsAttributes attributes;
138
139 const double step = snappedPoints.size() > 0 ? 100.0 / snappedPoints.size() : 1;
140 for ( int i = 0; i < snappedPoints.size(); i++ )
141 {
142 if ( feedback->isCanceled() )
143 {
144 break;
145 }
146
147 idxStart = graph->findVertex( snappedPoints.at( i ) );
148 origPoint = points.at( i ).toString();
149
150 QgsGraphAnalyzer::dijkstra( graph.get(), idxStart, 0, &tree, &costs );
151
152 QgsMultiPointXY areaPoints;
153 QgsMultiPolylineXY lines;
154 QSet< int > vertices;
155
156 for ( int j = 0; j < costs.size(); j++ )
157 {
158 inboundEdgeIndex = tree.at( j );
159
160 if ( inboundEdgeIndex == -1 && j != idxStart )
161 {
162 // unreachable vertex
163 continue;
164 }
165
166 startVertexCost = costs.at( j );
167 if ( startVertexCost > travelCost )
168 {
169 // vertex is too expensive, discard
170 continue;
171 }
172
173 vertices.insert( j );
174 startPoint = graph->vertex( j ).point();
175
176 // find all edges coming from this vertex
177 const QList< int > outgoingEdges = graph->vertex( j ).outgoingEdges() ;
178 for ( int edgeId : outgoingEdges )
179 {
180 edge = graph->edge( edgeId );
181 endVertexCost = startVertexCost + edge.cost( 0 ).toDouble();
182 endPoint = graph->vertex( edge.toVertex() ).point();
183 if ( endVertexCost <= travelCost )
184 {
185 // end vertex is cheap enough to include
186 vertices.insert( edge.toVertex() );
187 lines.push_back( QgsPolylineXY() << startPoint << endPoint );
188 }
189 else
190 {
191 // travelCost sits somewhere on this edge, interpolate position
192 QgsPointXY interpolatedEndPoint = QgsGeometryUtils::interpolatePointOnLineByValue( startPoint.x(), startPoint.y(), startVertexCost,
193 endPoint.x(), endPoint.y(), endVertexCost, travelCost );
194
195 areaPoints.push_back( interpolatedEndPoint );
196 lines.push_back( QgsPolylineXY() << startPoint << interpolatedEndPoint );
197 }
198 } // edges
199 } // costs
200
201 // convert to list and sort to maintain same order of points between algorithm runs
202 QList< int > verticesList = qgis::setToList( vertices );
203 areaPoints.reserve( verticesList.size() );
204 std::sort( verticesList.begin(), verticesList.end() );
205 for ( int v : verticesList )
206 {
207 areaPoints.push_back( graph->vertex( v ).point() );
208 }
209
210 if ( pointsSink )
211 {
212 QgsGeometry geomPoints = QgsGeometry::fromMultiPointXY( areaPoints );
213 feat.setGeometry( geomPoints );
214 attributes = sourceAttributes.value( i + 1 );
215 attributes << QStringLiteral( "within" ) << origPoint;
216 feat.setAttributes( attributes );
217 if ( !pointsSink->addFeature( feat, QgsFeatureSink::FastInsert ) )
218 throw QgsProcessingException( writeFeatureError( pointsSink.get(), parameters, QStringLiteral( "OUTPUT" ) ) );
219
220 if ( includeBounds )
221 {
222 QgsMultiPointXY upperBoundary, lowerBoundary;
223 QVector< int > nodes;
224 nodes.reserve( costs.size() );
225
226 int vertexId;
227 for ( int v = 0; v < costs.size(); v++ )
228 {
229 if ( costs.at( v ) > travelCost && tree.at( v ) != -1 )
230 {
231 vertexId = graph->edge( tree.at( v ) ).fromVertex();
232 if ( costs.at( vertexId ) <= travelCost )
233 {
234 nodes.push_back( v );
235 }
236 }
237 } // costs
238
239 upperBoundary.reserve( nodes.size() );
240 lowerBoundary.reserve( nodes.size() );
241 for ( int n : std::as_const( nodes ) )
242 {
243 upperBoundary.push_back( graph->vertex( graph->edge( tree.at( n ) ).toVertex() ).point() );
244 lowerBoundary.push_back( graph->vertex( graph->edge( tree.at( n ) ).fromVertex() ).point() );
245 } // nodes
246
247 QgsGeometry geomUpper = QgsGeometry::fromMultiPointXY( upperBoundary );
248 QgsGeometry geomLower = QgsGeometry::fromMultiPointXY( lowerBoundary );
249
250 feat.setGeometry( geomUpper );
251 attributes = sourceAttributes.value( i + 1 );
252 attributes << QStringLiteral( "upper" ) << origPoint;
253 feat.setAttributes( attributes );
254 if ( !pointsSink->addFeature( feat, QgsFeatureSink::FastInsert ) )
255 throw QgsProcessingException( writeFeatureError( pointsSink.get(), parameters, QStringLiteral( "OUTPUT" ) ) );
256
257 feat.setGeometry( geomLower );
258 attributes = sourceAttributes.value( i + 1 );
259 attributes << QStringLiteral( "lower" ) << origPoint;
260 feat.setAttributes( attributes );
261 if ( !pointsSink->addFeature( feat, QgsFeatureSink::FastInsert ) )
262 throw QgsProcessingException( writeFeatureError( pointsSink.get(), parameters, QStringLiteral( "OUTPUT" ) ) );
263 } // includeBounds
264 }
265
266 if ( linesSink )
267 {
269 feat.setGeometry( geomLines );
270 attributes = sourceAttributes.value( i + 1 );
271 attributes << QStringLiteral( "lines" ) << origPoint;
272 feat.setAttributes( attributes );
273 if ( !linesSink->addFeature( feat, QgsFeatureSink::FastInsert ) )
274 throw QgsProcessingException( writeFeatureError( linesSink.get(), parameters, QStringLiteral( "OUTPUT_LINES" ) ) );
275 }
276
277 feedback->setProgress( i * step );
278 } // snappedPoints
279
280 QVariantMap outputs;
281 if ( pointsSink )
282 {
283 outputs.insert( QStringLiteral( "OUTPUT" ), pointsSinkId );
284 }
285 if ( linesSink )
286 {
287 outputs.insert( QStringLiteral( "OUTPUT_LINES" ), linesSinkId );
288 }
289
290 return outputs;
291}
292
A vector of attributes.
Definition: qgsattributes.h:59
@ FastInsert
Use faster inserts, at the cost of updating the passed features to reflect changes made at the provid...
The feature class encapsulates a single feature including its unique ID, geometry and a list of field...
Definition: qgsfeature.h:56
void setAttributes(const QgsAttributes &attrs)
Sets the feature's attributes.
Definition: qgsfeature.cpp:160
void setGeometry(const QgsGeometry &geometry)
Set the feature's geometry.
Definition: qgsfeature.cpp:170
bool isCanceled() const SIP_HOLDGIL
Tells whether the operation has been canceled already.
Definition: qgsfeedback.h:54
void setProgress(double progress)
Sets the current progress for the feedback object.
Definition: qgsfeedback.h:63
Encapsulate a field in an attribute table or data source.
Definition: qgsfield.h:51
Container of fields for a vector layer.
Definition: qgsfields.h:45
bool append(const QgsField &field, FieldOrigin origin=OriginProvider, int originIndex=-1)
Appends a field. The field must have unique name, otherwise it is rejected (returns false)
Definition: qgsfields.cpp:59
static QgsPointXY interpolatePointOnLineByValue(double x1, double y1, double v1, double x2, double y2, double v2, double value) SIP_HOLDGIL
Interpolates the position of a point along the line from (x1, y1) to (x2, y2).
A geometry is the spatial representation of a feature.
Definition: qgsgeometry.h:164
static QgsGeometry fromMultiPolylineXY(const QgsMultiPolylineXY &multiline)
Creates a new geometry from a QgsMultiPolylineXY object.
static QgsGeometry fromMultiPointXY(const QgsMultiPointXY &multipoint)
Creates a new geometry from a QgsMultiPointXY object.
static void dijkstra(const QgsGraph *source, int startVertexIdx, int criterionNum, QVector< int > *resultTree=nullptr, QVector< double > *resultCost=nullptr)
Solve shortest path problem using Dijkstra algorithm.
This class implements a graph edge.
Definition: qgsgraph.h:45
int toVertex() const
Returns the index of the vertex at the end of this edge.
Definition: qgsgraph.cpp:184
QVariant cost(int strategyIndex) const
Returns edge cost calculated using specified strategy.
Definition: qgsgraph.cpp:169
A class to represent a 2D point.
Definition: qgspointxy.h:59
double y
Definition: qgspointxy.h:63
Q_GADGET double x
Definition: qgspointxy.h:62
Contains information about the context in which a processing algorithm is executed.
Custom exception class for processing related exceptions.
Definition: qgsexception.h:83
Base class for providing feedback from a processing algorithm.
virtual void pushInfo(const QString &info)
Pushes a general informational message from the algorithm.
@ FlagAdvanced
Parameter is an advanced parameter which should be hidden from users by default.
@ FlagHidden
Parameter is hidden and should not be shown to users.
An input feature source (such as vector layers) parameter for processing algorithms.
A numeric parameter for processing algorithms.
@ TypeVectorLine
Vector line layers.
Definition: qgsprocessing.h:50
@ TypeVectorPoint
Vector point layers.
Definition: qgsprocessing.h:49
QVector< QgsPolylineXY > QgsMultiPolylineXY
A collection of QgsPolylines that share a common collection of attributes.
Definition: qgsgeometry.h:86
QVector< QgsPointXY > QgsMultiPointXY
A collection of QgsPoints that share a common collection of attributes.
Definition: qgsgeometry.h:82
QVector< QgsPointXY > QgsPolylineXY
Polyline as represented as a vector of two-dimensional points.
Definition: qgsgeometry.h:63