QGIS API Documentation 3.38.0-Grenoble (exported)
Loading...
Searching...
No Matches
qgsalgorithmserviceareafrompoint.cpp
Go to the documentation of this file.
1/***************************************************************************
2 qgsalgorithmserviceareafrompoint.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 QgsServiceAreaFromPointAlgorithm::name() const
26{
27 return QStringLiteral( "serviceareafrompoint" );
28}
29
30QString QgsServiceAreaFromPointAlgorithm::displayName() const
31{
32 return QObject::tr( "Service area (from point)" );
33}
34
35QStringList QgsServiceAreaFromPointAlgorithm::tags() const
36{
37 return QObject::tr( "network,service,area,shortest,fastest" ).split( ',' );
38}
39
40QString QgsServiceAreaFromPointAlgorithm::shortHelpString() const
41{
42 return QObject::tr( "This algorithm creates a new vector with all the edges or parts of edges "
43 "of a network line layer that can be reached within a distance or a time, "
44 "starting from a point feature. The distance and the time (both referred to "
45 "as \"travel cost\") must be specified respectively in the network layer "
46 "units or in hours." );
47}
48
49QgsServiceAreaFromPointAlgorithm *QgsServiceAreaFromPointAlgorithm::createInstance() const
50{
51 return new QgsServiceAreaFromPointAlgorithm();
52}
53
54void QgsServiceAreaFromPointAlgorithm::initAlgorithm( const QVariantMap & )
55{
56 addCommonParams();
57 addParameter( new QgsProcessingParameterPoint( QStringLiteral( "START_POINT" ), QObject::tr( "Start point" ) ) );
58
59 std::unique_ptr< QgsProcessingParameterNumber > travelCost = std::make_unique< QgsProcessingParameterNumber >( QStringLiteral( "TRAVEL_COST" ), QObject::tr( "Travel cost (distance for 'Shortest', time for 'Fastest')" ), Qgis::ProcessingNumberParameterType::Double, 0, true, 0 );
60 travelCost->setFlags( travelCost->flags() | Qgis::ProcessingParameterFlag::Hidden );
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< QgsProcessingParameterNumber > maxPointDistanceFromNetwork = std::make_unique < QgsProcessingParameterDistance >( QStringLiteral( "POINT_TOLERANCE" ), QObject::tr( "Maximum point distance from network" ), QVariant(), QStringLiteral( "INPUT" ), true, 0 );
67 maxPointDistanceFromNetwork->setFlags( maxPointDistanceFromNetwork->flags() | Qgis::ProcessingParameterFlag::Advanced );
68 maxPointDistanceFromNetwork->setHelp( QObject::tr( "Specifies an optional limit on the distance from the point to the network layer. If the point is further from the network than this distance an error will be raised." ) );
69 addParameter( maxPointDistanceFromNetwork.release() );
70
71 std::unique_ptr< QgsProcessingParameterBoolean > includeBounds = std::make_unique< QgsProcessingParameterBoolean >( QStringLiteral( "INCLUDE_BOUNDS" ), QObject::tr( "Include upper/lower bound points" ), false, true );
72 includeBounds->setFlags( includeBounds->flags() | Qgis::ProcessingParameterFlag::Advanced );
73 addParameter( includeBounds.release() );
74
75 std::unique_ptr< QgsProcessingParameterFeatureSink > outputLines = std::make_unique< QgsProcessingParameterFeatureSink >( QStringLiteral( "OUTPUT_LINES" ), QObject::tr( "Service area (lines)" ),
77 outputLines->setCreateByDefault( true );
78 addParameter( outputLines.release() );
79
80 std::unique_ptr< QgsProcessingParameterFeatureSink > outputPoints = std::make_unique< QgsProcessingParameterFeatureSink >( QStringLiteral( "OUTPUT" ), QObject::tr( "Service area (boundary nodes)" ),
82 outputPoints->setCreateByDefault( false );
83 addParameter( outputPoints.release() );
84}
85
86QVariantMap QgsServiceAreaFromPointAlgorithm::processAlgorithm( const QVariantMap &parameters, QgsProcessingContext &context, QgsProcessingFeedback *feedback )
87{
88 loadCommonParams( parameters, context, feedback );
89
90 const QgsPointXY startPoint = parameterAsPoint( parameters, QStringLiteral( "START_POINT" ), context, mNetwork->sourceCrs() );
91
92 // use older deprecated travel cost style if specified, to maintain old api
93 const bool useOldTravelCost = parameters.value( QStringLiteral( "TRAVEL_COST" ) ).isValid();
94 double travelCost = parameterAsDouble( parameters, useOldTravelCost ? QStringLiteral( "TRAVEL_COST" ) : QStringLiteral( "TRAVEL_COST2" ), context );
95
96 const int strategy = parameterAsInt( parameters, QStringLiteral( "STRATEGY" ), context );
97 if ( strategy && !useOldTravelCost )
98 travelCost *= mMultiplier;
99
100 bool includeBounds = true; // default to true to maintain 3.0 API
101 if ( parameters.contains( QStringLiteral( "INCLUDE_BOUNDS" ) ) )
102 {
103 includeBounds = parameterAsBool( parameters, QStringLiteral( "INCLUDE_BOUNDS" ), context );
104 }
105
106 feedback->pushInfo( QObject::tr( "Building graph…" ) );
107 QVector< QgsPointXY > snappedPoints;
108 mDirector->makeGraph( mBuilder.get(), { startPoint }, snappedPoints, feedback );
109 const QgsPointXY snappedStartPoint = snappedPoints[0];
110
111 // check distance for the snapped point
112 if ( parameters.value( QStringLiteral( "POINT_TOLERANCE" ) ).isValid() )
113 {
114 const double pointDistanceThreshold = parameterAsDouble( parameters, QStringLiteral( "POINT_TOLERANCE" ), context );
115
116 const double distancePointToNetwork = mBuilder->distanceArea()->measureLine( startPoint, snappedStartPoint );
117 if ( distancePointToNetwork > pointDistanceThreshold )
118 {
119 throw QgsProcessingException( QObject::tr( "Point is too far from the network layer (%1, maximum permitted is %2)" ).arg( distancePointToNetwork ).arg( pointDistanceThreshold ) );
120 }
121 }
122
123 feedback->pushInfo( QObject::tr( "Calculating service area…" ) );
124 std::unique_ptr< QgsGraph> graph( mBuilder->takeGraph() );
125 const int idxStart = graph->findVertex( snappedStartPoint );
126
127 QVector< int > tree;
128 QVector< double > costs;
129 QgsGraphAnalyzer::dijkstra( graph.get(), idxStart, 0, &tree, &costs );
130
131 QgsMultiPointXY points;
132 QgsMultiPolylineXY lines;
133 QSet< int > vertices;
134
135 int inboundEdgeIndex;
136 double startVertexCost, endVertexCost;
137 QgsPointXY edgeStart, edgeEnd;
138 QgsGraphEdge edge;
139
140 for ( int i = 0; i < costs.size(); i++ )
141 {
142 inboundEdgeIndex = tree.at( i );
143 if ( inboundEdgeIndex == -1 && i != idxStart )
144 {
145 // unreachable vertex
146 continue;
147 }
148
149 startVertexCost = costs.at( i );
150 if ( startVertexCost > travelCost )
151 {
152 // vertex is too expensive, discard
153 continue;
154 }
155
156 vertices.insert( i );
157 edgeStart = graph->vertex( i ).point();
158
159 // find all edges coming from this vertex
160 const QList< int > outgoingEdges = graph->vertex( i ).outgoingEdges() ;
161 for ( const int edgeId : outgoingEdges )
162 {
163 edge = graph->edge( edgeId );
164 endVertexCost = startVertexCost + edge.cost( 0 ).toDouble();
165 edgeEnd = graph->vertex( edge.toVertex() ).point();
166 if ( endVertexCost <= travelCost )
167 {
168 // end vertex is cheap enough to include
169 vertices.insert( edge.toVertex() );
170 lines.push_back( QgsPolylineXY() << edgeStart << edgeEnd );
171 }
172 else
173 {
174 // travelCost sits somewhere on this edge, interpolate position
175 const QgsPointXY interpolatedEndPoint = QgsGeometryUtils::interpolatePointOnLineByValue( edgeStart.x(), edgeStart.y(), startVertexCost,
176 edgeEnd.x(), edgeEnd.y(), endVertexCost, travelCost );
177
178 points.push_back( interpolatedEndPoint );
179 lines.push_back( QgsPolylineXY() << edgeStart << interpolatedEndPoint );
180 }
181 } // edges
182 } // costs
183
184 // convert to list and sort to maintain same order of points between algorithm runs
185 QList< int > verticesList = qgis::setToList( vertices );
186 points.reserve( verticesList.size() );
187 std::sort( verticesList.begin(), verticesList.end() );
188 for ( const int v : verticesList )
189 {
190 points.push_back( graph->vertex( v ).point() );
191 }
192
193 feedback->pushInfo( QObject::tr( "Writing results…" ) );
194
195 QVariantMap outputs;
196
197 QgsFields fields;
198 fields.append( QgsField( QStringLiteral( "type" ), QMetaType::Type::QString ) );
199 fields.append( QgsField( QStringLiteral( "start" ), QMetaType::Type::QString ) );
200
201 QgsFeature feat;
202 feat.setFields( fields );
203
204 QString pointsSinkId;
205 std::unique_ptr< QgsFeatureSink > pointsSink( parameterAsSink( parameters, QStringLiteral( "OUTPUT" ), context, pointsSinkId, fields,
206 Qgis::WkbType::MultiPoint, mNetwork->sourceCrs() ) );
207
208 if ( pointsSink )
209 {
210 outputs.insert( QStringLiteral( "OUTPUT" ), pointsSinkId );
211
212 const QgsGeometry geomPoints = QgsGeometry::fromMultiPointXY( points );
213 feat.setGeometry( geomPoints );
214 feat.setAttributes( QgsAttributes() << QStringLiteral( "within" ) << startPoint.toString() );
215 if ( !pointsSink->addFeature( feat, QgsFeatureSink::FastInsert ) )
216 throw QgsProcessingException( writeFeatureError( pointsSink.get(), parameters, QStringLiteral( "OUTPUT" ) ) );
217
218 if ( includeBounds )
219 {
220 QgsMultiPointXY upperBoundary, lowerBoundary;
221 QVector< int > nodes;
222
223 int vertexId;
224 for ( int i = 0; i < costs.size(); i++ )
225 {
226 if ( costs.at( i ) > travelCost && tree.at( i ) != -1 )
227 {
228 vertexId = graph->edge( tree.at( i ) ).fromVertex();
229 if ( costs.at( vertexId ) <= travelCost )
230 {
231 nodes.push_back( i );
232 }
233 }
234 } // costs
235
236 upperBoundary.reserve( nodes.size() );
237 lowerBoundary.reserve( nodes.size() );
238 for ( const int i : nodes )
239 {
240 upperBoundary.push_back( graph->vertex( graph->edge( tree.at( i ) ).toVertex() ).point() );
241 lowerBoundary.push_back( graph->vertex( graph->edge( tree.at( i ) ).fromVertex() ).point() );
242 } // nodes
243
244 const QgsGeometry geomUpper = QgsGeometry::fromMultiPointXY( upperBoundary );
245 const QgsGeometry geomLower = QgsGeometry::fromMultiPointXY( lowerBoundary );
246
247 feat.setGeometry( geomUpper );
248 feat.setAttributes( QgsAttributes() << QStringLiteral( "upper" ) << startPoint.toString() );
249 if ( !pointsSink->addFeature( feat, QgsFeatureSink::FastInsert ) )
250 throw QgsProcessingException( writeFeatureError( pointsSink.get(), parameters, QStringLiteral( "OUTPUT" ) ) );
251
252 feat.setGeometry( geomLower );
253 feat.setAttributes( QgsAttributes() << QStringLiteral( "lower" ) << startPoint.toString() );
254 if ( !pointsSink->addFeature( feat, QgsFeatureSink::FastInsert ) )
255 throw QgsProcessingException( writeFeatureError( pointsSink.get(), parameters, QStringLiteral( "OUTPUT" ) ) );
256 } // includeBounds
257 }
258
259 QString linesSinkId;
260 std::unique_ptr< QgsFeatureSink > linesSink( parameterAsSink( parameters, QStringLiteral( "OUTPUT_LINES" ), context, linesSinkId, fields,
261 Qgis::WkbType::MultiLineString, mNetwork->sourceCrs() ) );
262
263 if ( linesSink )
264 {
265 outputs.insert( QStringLiteral( "OUTPUT_LINES" ), linesSinkId );
266 const QgsGeometry geomLines = QgsGeometry::fromMultiPolylineXY( lines );
267 feat.setGeometry( geomLines );
268 feat.setAttributes( QgsAttributes() << QStringLiteral( "lines" ) << startPoint.toString() );
269 if ( !linesSink->addFeature( feat, QgsFeatureSink::FastInsert ) )
270 throw QgsProcessingException( writeFeatureError( linesSink.get(), parameters, QStringLiteral( "OUTPUT_LINES" ) ) );
271 }
272
273 return outputs;
274}
275
@ VectorPoint
Vector point layers.
@ VectorLine
Vector line layers.
@ MultiPoint
MultiPoint.
@ MultiLineString
MultiLineString.
@ Hidden
Parameter is hidden and should not be shown to users.
@ Advanced
Parameter is an advanced parameter which should be hidden from users by default.
A vector of attributes.
@ 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:58
void setAttributes(const QgsAttributes &attrs)
Sets the feature's attributes.
void setFields(const QgsFields &fields, bool initAttributes=false)
Assigns a field map with the feature to allow attribute access by attribute name.
void setGeometry(const QgsGeometry &geometry)
Set the feature's geometry.
Encapsulate a field in an attribute table or data source.
Definition qgsfield.h:53
Container of fields for a vector layer.
Definition qgsfields.h:46
bool append(const QgsField &field, Qgis::FieldOrigin origin=Qgis::FieldOrigin::Provider, int originIndex=-1)
Appends a field.
Definition qgsfields.cpp:60
static QgsPointXY interpolatePointOnLineByValue(double x1, double y1, double v1, double x2, double y2, double v2, double value)
Interpolates the position of a point along the line from (x1, y1) to (x2, y2).
A geometry is the spatial representation of a feature.
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:44
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:60
QString toString(int precision=-1) const
Returns a string representation of the point (x, y) with a preset precision.
double y
Definition qgspointxy.h:64
double x
Definition qgspointxy.h:63
Contains information about the context in which a processing algorithm is executed.
Custom exception class for processing related exceptions.
Base class for providing feedback from a processing algorithm.
virtual void pushInfo(const QString &info)
Pushes a general informational message from the algorithm.
A numeric parameter for processing algorithms.
A point parameter for processing algorithms.
QVector< QgsPolylineXY > QgsMultiPolylineXY
A collection of QgsPolylines that share a common collection of attributes.
Definition qgsgeometry.h:84
QVector< QgsPointXY > QgsMultiPointXY
A collection of QgsPoints that share a common collection of attributes.
Definition qgsgeometry.h:80
QVector< QgsPointXY > QgsPolylineXY
Polyline as represented as a vector of two-dimensional points.
Definition qgsgeometry.h:62