QGIS API Documentation 3.40.0-Bratislava (b56115d8743)
Loading...
Searching...
No Matches
qgsalgorithmjoinwithlines.cpp
Go to the documentation of this file.
1/***************************************************************************
2 qgsalgorithmjoinwithlines.cpp
3 ---------------------
4 begin : April 2017
5 copyright : (C) 2017 by Nyall Dawson
6 email : nyall dot dawson 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#include "qgslinestring.h"
20#include "qgsmultilinestring.h"
21#include "qgsdistancearea.h"
22
24
25QString QgsJoinWithLinesAlgorithm::name() const
26{
27 return QStringLiteral( "hublines" );
28}
29
30QString QgsJoinWithLinesAlgorithm::displayName() const
31{
32 return QObject::tr( "Join by lines (hub lines)" );
33}
34
35QStringList QgsJoinWithLinesAlgorithm::tags() const
36{
37 return QObject::tr( "join,connect,lines,points,hub,spoke,geodesic,great,circle" ).split( ',' );
38}
39
40QString QgsJoinWithLinesAlgorithm::group() const
41{
42 return QObject::tr( "Vector analysis" );
43}
44
45QString QgsJoinWithLinesAlgorithm::groupId() const
46{
47 return QStringLiteral( "vectoranalysis" );
48}
49
50void QgsJoinWithLinesAlgorithm::initAlgorithm( const QVariantMap & )
51{
52 addParameter( new QgsProcessingParameterFeatureSource( QStringLiteral( "HUBS" ),
53 QObject::tr( "Hub layer" ) ) );
54 addParameter( new QgsProcessingParameterField( QStringLiteral( "HUB_FIELD" ),
55 QObject::tr( "Hub ID field" ), QVariant(), QStringLiteral( "HUBS" ) ) );
56
57 addParameter( new QgsProcessingParameterField( QStringLiteral( "HUB_FIELDS" ),
58 QObject::tr( "Hub layer fields to copy (leave empty to copy all fields)" ),
59 QVariant(), QStringLiteral( "HUBS" ), Qgis::ProcessingFieldParameterDataType::Any,
60 true, true ) );
61
62 addParameter( new QgsProcessingParameterFeatureSource( QStringLiteral( "SPOKES" ),
63 QObject::tr( "Spoke layer" ) ) );
64 addParameter( new QgsProcessingParameterField( QStringLiteral( "SPOKE_FIELD" ),
65 QObject::tr( "Spoke ID field" ), QVariant(), QStringLiteral( "SPOKES" ) ) );
66
67 addParameter( new QgsProcessingParameterField( QStringLiteral( "SPOKE_FIELDS" ),
68 QObject::tr( "Spoke layer fields to copy (leave empty to copy all fields)" ),
69 QVariant(), QStringLiteral( "SPOKES" ), Qgis::ProcessingFieldParameterDataType::Any,
70 true, true ) );
71
72 addParameter( new QgsProcessingParameterBoolean( QStringLiteral( "GEODESIC" ), QObject::tr( "Create geodesic lines" ), false ) );
73
74 auto distanceParam = std::make_unique< QgsProcessingParameterDistance >( QStringLiteral( "GEODESIC_DISTANCE" ), QObject::tr( "Distance between vertices (geodesic lines only)" ), 1000 );
75 distanceParam->setFlags( distanceParam->flags() | Qgis::ProcessingParameterFlag::Advanced );
76 distanceParam->setDefaultUnit( Qgis::DistanceUnit::Kilometers );
77 distanceParam->setIsDynamic( true );
78 distanceParam->setDynamicPropertyDefinition( QgsPropertyDefinition( QStringLiteral( "Geodesic Distance" ), QObject::tr( "Distance between vertices" ), QgsPropertyDefinition::DoublePositive ) );
79 distanceParam->setDynamicLayerParameterName( QStringLiteral( "HUBS" ) );
80 addParameter( distanceParam.release() );
81
82 auto breakParam = std::make_unique< QgsProcessingParameterBoolean >( QStringLiteral( "ANTIMERIDIAN_SPLIT" ), QObject::tr( "Split lines at antimeridian (±180 degrees longitude)" ), false );
83 breakParam->setFlags( breakParam->flags() | Qgis::ProcessingParameterFlag::Advanced );
84 addParameter( breakParam.release() );
85
86 addParameter( new QgsProcessingParameterFeatureSink( QStringLiteral( "OUTPUT" ), QObject::tr( "Hub lines" ), Qgis::ProcessingSourceType::VectorLine ) );
87}
88
89QString QgsJoinWithLinesAlgorithm::shortHelpString() const
90{
91 return QObject::tr( "This algorithm creates hub and spoke diagrams by connecting lines from points on the Spoke layer to matching points in the Hub layer.\n\n"
92 "Determination of which hub goes with each point is based on a match between the Hub ID field on the hub points and the Spoke ID field on the spoke points.\n\n"
93 "If input layers are not point layers, a point on the surface of the geometries will be taken as the connecting location.\n\n"
94 "Optionally, geodesic lines can be created, which represent the shortest path on the surface of an ellipsoid. When "
95 "geodesic mode is used, it is possible to split the created lines at the antimeridian (±180 degrees longitude), which can improve "
96 "rendering of the lines. Additionally, the distance between vertices can be specified. A smaller distance results in a denser, more "
97 "accurate line." );
98}
99
100QString QgsJoinWithLinesAlgorithm::shortDescription() const
101{
102 return QObject::tr( "Creates lines joining two point layers, based on a common attribute value." );
103}
104
105Qgis::ProcessingAlgorithmDocumentationFlags QgsJoinWithLinesAlgorithm::documentationFlags() const
106{
108}
109
110QgsJoinWithLinesAlgorithm *QgsJoinWithLinesAlgorithm::createInstance() const
111{
112 return new QgsJoinWithLinesAlgorithm();
113}
114
115QVariantMap QgsJoinWithLinesAlgorithm::processAlgorithm( const QVariantMap &parameters, QgsProcessingContext &context, QgsProcessingFeedback *feedback )
116{
117 if ( parameters.value( QStringLiteral( "SPOKES" ) ) == parameters.value( QStringLiteral( "HUBS" ) ) )
118 throw QgsProcessingException( QObject::tr( "Same layer given for both hubs and spokes" ) );
119
120 std::unique_ptr< QgsProcessingFeatureSource > hubSource( parameterAsSource( parameters, QStringLiteral( "HUBS" ), context ) );
121 if ( !hubSource )
122 throw QgsProcessingException( invalidSourceError( parameters, QStringLiteral( "HUBS" ) ) );
123
124 std::unique_ptr< QgsProcessingFeatureSource > spokeSource( parameterAsSource( parameters, QStringLiteral( "SPOKES" ), context ) );
125 if ( !spokeSource )
126 throw QgsProcessingException( invalidSourceError( parameters, QStringLiteral( "SPOKES" ) ) );
127
128 const QString fieldHubName = parameterAsString( parameters, QStringLiteral( "HUB_FIELD" ), context );
129 const int fieldHubIndex = hubSource->fields().lookupField( fieldHubName );
130 const QStringList hubFieldsToCopy = parameterAsStrings( parameters, QStringLiteral( "HUB_FIELDS" ), context );
131
132 const QString fieldSpokeName = parameterAsString( parameters, QStringLiteral( "SPOKE_FIELD" ), context );
133 const int fieldSpokeIndex = spokeSource->fields().lookupField( fieldSpokeName );
134 const QStringList spokeFieldsToCopy = parameterAsStrings( parameters, QStringLiteral( "SPOKE_FIELDS" ), context );
135
136 if ( fieldHubIndex < 0 || fieldSpokeIndex < 0 )
137 throw QgsProcessingException( QObject::tr( "Invalid ID field" ) );
138
139 const bool geodesic = parameterAsBoolean( parameters, QStringLiteral( "GEODESIC" ), context );
140 const double geodesicDistance = parameterAsDouble( parameters, QStringLiteral( "GEODESIC_DISTANCE" ), context ) * 1000;
141 const bool dynamicGeodesicDistance = QgsProcessingParameters::isDynamic( parameters, QStringLiteral( "GEODESIC_DISTANCE" ) );
142 QgsExpressionContext expressionContext = createExpressionContext( parameters, context, hubSource.get() );
143 QgsProperty geodesicDistanceProperty;
144 if ( dynamicGeodesicDistance )
145 {
146 geodesicDistanceProperty = parameters.value( QStringLiteral( "GEODESIC_DISTANCE" ) ).value< QgsProperty >();
147 }
148
149 const bool splitAntimeridian = parameterAsBoolean( parameters, QStringLiteral( "ANTIMERIDIAN_SPLIT" ), context );
151 da.setSourceCrs( hubSource->sourceCrs(), context.transformContext() );
152 da.setEllipsoid( context.ellipsoid() );
153
154 QgsFields hubOutFields;
155 QgsAttributeList hubFieldIndices;
156 if ( hubFieldsToCopy.empty() )
157 {
158 hubOutFields = hubSource->fields();
159 hubFieldIndices.reserve( hubOutFields.count() );
160 for ( int i = 0; i < hubOutFields.count(); ++i )
161 {
162 hubFieldIndices << i;
163 }
164 }
165 else
166 {
167 hubFieldIndices.reserve( hubOutFields.count() );
168 for ( const QString &field : hubFieldsToCopy )
169 {
170 const int index = hubSource->fields().lookupField( field );
171 if ( index >= 0 )
172 {
173 hubFieldIndices << index;
174 hubOutFields.append( hubSource->fields().at( index ) );
175 }
176 }
177 }
178
179 QgsAttributeList hubFields2Fetch = hubFieldIndices;
180 hubFields2Fetch << fieldHubIndex;
181
182 QgsFields spokeOutFields;
183 QgsAttributeList spokeFieldIndices;
184 if ( spokeFieldsToCopy.empty() )
185 {
186 spokeOutFields = spokeSource->fields();
187 spokeFieldIndices.reserve( spokeOutFields.count() );
188 for ( int i = 0; i < spokeOutFields.count(); ++i )
189 {
190 spokeFieldIndices << i;
191 }
192 }
193 else
194 {
195 for ( const QString &field : spokeFieldsToCopy )
196 {
197 const int index = spokeSource->fields().lookupField( field );
198 if ( index >= 0 )
199 {
200 spokeFieldIndices << index;
201 spokeOutFields.append( spokeSource->fields().at( index ) );
202 }
203 }
204 }
205
206 QgsAttributeList spokeFields2Fetch = spokeFieldIndices;
207 spokeFields2Fetch << fieldSpokeIndex;
208
209
210 const QgsFields fields = QgsProcessingUtils::combineFields( hubOutFields, spokeOutFields );
211
213 bool hasZ = false;
214 if ( !geodesic && ( QgsWkbTypes::hasZ( hubSource->wkbType() ) || QgsWkbTypes::hasZ( spokeSource->wkbType() ) ) )
215 {
216 outType = QgsWkbTypes::addZ( outType );
217 hasZ = true;
218 }
219 bool hasM = false;
220 if ( !geodesic && ( QgsWkbTypes::hasM( hubSource->wkbType() ) || QgsWkbTypes::hasM( spokeSource->wkbType() ) ) )
221 {
222 outType = QgsWkbTypes::addM( outType );
223 hasM = true;
224 }
225
226 QString dest;
227 std::unique_ptr< QgsFeatureSink > sink( parameterAsSink( parameters, QStringLiteral( "OUTPUT" ), context, dest, fields,
228 outType, hubSource->sourceCrs(), QgsFeatureSink::RegeneratePrimaryKey ) );
229 if ( !sink )
230 throw QgsProcessingException( invalidSinkError( parameters, QStringLiteral( "OUTPUT" ) ) );
231
232 auto getPointFromFeature = [hasZ, hasM]( const QgsFeature & feature )->QgsPoint
233 {
234 QgsPoint p;
235 if ( feature.geometry().type() == Qgis::GeometryType::Point && !feature.geometry().isMultipart() )
236 p = *static_cast< const QgsPoint *>( feature.geometry().constGet() );
237 else
238 p = *static_cast< const QgsPoint *>( feature.geometry().pointOnSurface().constGet() );
239 if ( hasZ && !p.is3D() )
240 p.addZValue( 0 );
241 if ( hasM && !p.isMeasure() )
242 p.addMValue( 0 );
243 return p;
244 };
245
247 const double step = hubSource->featureCount() > 0 ? 100.0 / hubSource->featureCount() : 1;
248 int i = 0;
249 QgsFeature hubFeature;
250 while ( hubFeatures.nextFeature( hubFeature ) )
251 {
252 i++;
253 if ( feedback->isCanceled() )
254 {
255 break;
256 }
257
258 feedback->setProgress( i * step );
259
260 if ( !hubFeature.hasGeometry() )
261 continue;
262
263 const QgsPoint hubPoint = getPointFromFeature( hubFeature );
264
265 // only keep selected attributes
266 QgsAttributes hubAttributes;
267 const int attributeCount = hubFeature.attributeCount();
268 for ( int j = 0; j < attributeCount; ++j )
269 {
270 if ( !hubFieldIndices.contains( j ) )
271 continue;
272 hubAttributes << hubFeature.attribute( j );
273 }
274
275 QgsFeatureRequest spokeRequest = QgsFeatureRequest().setDestinationCrs( hubSource->sourceCrs(), context.transformContext() );
276 spokeRequest.setSubsetOfAttributes( spokeFields2Fetch );
277 spokeRequest.setFilterExpression( QgsExpression::createFieldEqualityExpression( fieldSpokeName, hubFeature.attribute( fieldHubIndex ) ) );
278
279 QgsFeatureIterator spokeFeatures = spokeSource->getFeatures( spokeRequest, Qgis::ProcessingFeatureSourceFlag::SkipGeometryValidityChecks );
280 QgsFeature spokeFeature;
281 while ( spokeFeatures.nextFeature( spokeFeature ) )
282 {
283 if ( feedback->isCanceled() )
284 {
285 break;
286 }
287 if ( !spokeFeature.hasGeometry() )
288 continue;
289
290 const QgsPoint spokePoint = getPointFromFeature( spokeFeature );
291 QgsGeometry line;
292 if ( !geodesic )
293 {
294 line = QgsGeometry( new QgsLineString( QVector< QgsPoint >() << hubPoint << spokePoint ) );
295 if ( splitAntimeridian )
296 line = da.splitGeometryAtAntimeridian( line );
297 }
298 else
299 {
300 double distance = geodesicDistance;
301 if ( dynamicGeodesicDistance )
302 {
303 expressionContext.setFeature( hubFeature );
304 distance = geodesicDistanceProperty.valueAsDouble( expressionContext, distance );
305 }
306
307 std::unique_ptr< QgsMultiLineString > ml = std::make_unique< QgsMultiLineString >();
308 std::unique_ptr< QgsLineString > l = std::make_unique< QgsLineString >( QVector< QgsPoint >() << hubPoint );
309 const QVector< QVector< QgsPointXY > > points = da.geodesicLine( QgsPointXY( hubPoint ), QgsPointXY( spokePoint ), distance, splitAntimeridian );
310 QVector< QgsPointXY > points1 = points.at( 0 );
311 points1.pop_front();
312 if ( points.count() == 1 )
313 points1.pop_back();
314
315 const QgsLineString geodesicPoints( points1 );
316 l->append( &geodesicPoints );
317 if ( points.count() == 1 )
318 l->addVertex( spokePoint );
319
320 ml->addGeometry( l.release() );
321 if ( points.count() > 1 )
322 {
323 QVector< QgsPointXY > points2 = points.at( 1 );
324 points2.pop_back();
325 l = std::make_unique< QgsLineString >( points2 );
326 if ( hasZ )
327 l->addZValue( std::numeric_limits<double>::quiet_NaN() );
328 if ( hasM )
329 l->addMValue( std::numeric_limits<double>::quiet_NaN() );
330
331 l->addVertex( spokePoint );
332 ml->addGeometry( l.release() );
333 }
334 line = QgsGeometry( std::move( ml ) );
335 }
336
337 QgsFeature outFeature;
338 QgsAttributes outAttributes = hubAttributes;
339
340 // only keep selected attributes
341 QgsAttributes spokeAttributes;
342 const int attributeCount = spokeFeature.attributeCount();
343 for ( int j = 0; j < attributeCount; ++j )
344 {
345 if ( !spokeFieldIndices.contains( j ) )
346 continue;
347 spokeAttributes << spokeFeature.attribute( j );
348 }
349
350 outAttributes.append( spokeAttributes );
351 outFeature.setAttributes( outAttributes );
352 outFeature.setGeometry( line );
353 if ( !sink->addFeature( outFeature, QgsFeatureSink::FastInsert ) )
354 throw QgsProcessingException( writeFeatureError( sink.get(), parameters, QStringLiteral( "OUTPUT" ) ) );
355 }
356 }
357
358 QVariantMap outputs;
359 outputs.insert( QStringLiteral( "OUTPUT" ), dest );
360 return outputs;
361}
362
@ VectorLine
Vector line layers.
@ Kilometers
Kilometers.
@ RegeneratesPrimaryKey
Algorithm always drops any existing primary keys or FID values and regenerates them in outputs.
QFlags< ProcessingAlgorithmDocumentationFlag > ProcessingAlgorithmDocumentationFlags
Flags describing algorithm behavior for documentation purposes.
Definition qgis.h:3359
@ SkipGeometryValidityChecks
Invalid geometry checks should always be skipped. This flag can be useful for algorithms which always...
WkbType
The WKB type describes the number of dimensions a geometry has.
Definition qgis.h:256
@ LineString
LineString.
@ MultiLineString
MultiLineString.
@ Advanced
Parameter is an advanced parameter which should be hidden from users by default.
bool isMeasure() const
Returns true if the geometry contains m values.
bool is3D() const
Returns true if the geometry is 3D and contains a z-value.
A vector of attributes.
A general purpose distance and area calculator, capable of performing ellipsoid based calculations.
QVector< QVector< QgsPointXY > > geodesicLine(const QgsPointXY &p1, const QgsPointXY &p2, double interval, bool breakLine=false) const
Calculates the geodesic line between p1 and p2, which represents the shortest path on the ellipsoid b...
void setSourceCrs(const QgsCoordinateReferenceSystem &crs, const QgsCoordinateTransformContext &context)
Sets source spatial reference system crs.
QgsGeometry splitGeometryAtAntimeridian(const QgsGeometry &geometry) const
Splits a (Multi)LineString geometry at the antimeridian (longitude +/- 180 degrees).
bool setEllipsoid(const QString &ellipsoid)
Sets the ellipsoid by its acronym.
Expression contexts are used to encapsulate the parameters around which a QgsExpression should be eva...
void setFeature(const QgsFeature &feature)
Convenience function for setting a feature for the context.
static QString createFieldEqualityExpression(const QString &fieldName, const QVariant &value, QMetaType::Type fieldType=QMetaType::Type::UnknownType)
Create an expression allowing to evaluate if a field is equal to a value.
Wrapper for iterator of features from vector data provider or vector layer.
bool nextFeature(QgsFeature &f)
Fetch next feature and stores in f, returns true on success.
This class wraps a request for features to a vector layer (or directly its vector data provider).
QgsFeatureRequest & setSubsetOfAttributes(const QgsAttributeList &attrs)
Set a subset of attributes that will be fetched.
QgsFeatureRequest & setDestinationCrs(const QgsCoordinateReferenceSystem &crs, const QgsCoordinateTransformContext &context)
Sets the destination crs for feature's geometries.
QgsFeatureRequest & setFilterExpression(const QString &expression)
Set the filter expression.
@ FastInsert
Use faster inserts, at the cost of updating the passed features to reflect changes made at the provid...
@ RegeneratePrimaryKey
This flag indicates, that a primary key field cannot be guaranteed to be unique and the sink should i...
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.
int attributeCount() const
Returns the number of attributes attached to the feature.
bool hasGeometry() const
Returns true if the feature has an associated geometry.
Q_INVOKABLE QVariant attribute(const QString &name) const
Lookup attribute value by attribute name.
void setGeometry(const QgsGeometry &geometry)
Set the feature's geometry.
bool isCanceled() const
Tells whether the operation has been canceled already.
Definition qgsfeedback.h:53
void setProgress(double progress)
Sets the current progress for the feedback object.
Definition qgsfeedback.h:61
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:69
int count
Definition qgsfields.h:50
Q_INVOKABLE int lookupField(const QString &fieldName) const
Looks up field's index from the field name.
A geometry is the spatial representation of a feature.
Line string geometry type, with support for z-dimension and m-values.
A class to represent a 2D point.
Definition qgspointxy.h:60
Point geometry type, with support for z-dimension and m-values.
Definition qgspoint.h:49
bool addMValue(double mValue=0) override
Adds a measure to the geometry, initialized to a preset value.
Definition qgspoint.cpp:568
bool addZValue(double zValue=0) override
Adds a z-dimension to the geometry, initialized to a preset value.
Definition qgspoint.cpp:557
Contains information about the context in which a processing algorithm is executed.
QgsCoordinateTransformContext transformContext() const
Returns the coordinate transform context.
QString ellipsoid() const
Returns the ellipsoid to use for distance and area calculations.
Custom exception class for processing related exceptions.
Base class for providing feedback from a processing algorithm.
A boolean parameter for processing algorithms.
A feature sink output for processing algorithms.
An input feature source (such as vector layers) parameter for processing algorithms.
A vector layer or feature source field parameter for processing algorithms.
static bool isDynamic(const QVariantMap &parameters, const QString &name)
Returns true if the parameter with matching name is a dynamic parameter, and must be evaluated once f...
static QgsFields combineFields(const QgsFields &fieldsA, const QgsFields &fieldsB, const QString &fieldsBPrefix=QString())
Combines two field lists, avoiding duplicate field names (in a case-insensitive manner).
Definition for a property.
Definition qgsproperty.h:45
@ DoublePositive
Positive double value (including 0)
Definition qgsproperty.h:56
A store for object properties.
QVariant value(const QgsExpressionContext &context, const QVariant &defaultValue=QVariant(), bool *ok=nullptr) const
Calculates the current value of the property, including any transforms which are set for the property...
double valueAsDouble(const QgsExpressionContext &context, double defaultValue=0.0, bool *ok=nullptr) const
Calculates the current value of the property and interprets it as a double.
static Qgis::WkbType addM(Qgis::WkbType type)
Adds the m dimension to a WKB type and returns the new type.
static Qgis::WkbType addZ(Qgis::WkbType type)
Adds the z dimension to a WKB type and returns the new type.
static bool hasZ(Qgis::WkbType type)
Tests whether a WKB type contains the z-dimension.
static bool hasM(Qgis::WkbType type)
Tests whether a WKB type contains m values.
QList< int > QgsAttributeList
Definition qgsfield.h:27