QGIS API Documentation 3.99.0-Master (2fe06baccd8)
Loading...
Searching...
No Matches
qgsalgorithmexportmesh.cpp
Go to the documentation of this file.
1/***************************************************************************
2 qgsalgorithmexportmesh.cpp
3 ---------------------------
4 begin : October 2020
5 copyright : (C) 2020 by Vincent Cloarec
6 email : vcloarec 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 "qgslinestring.h"
21#include "qgsmeshcontours.h"
22#include "qgsmeshdataset.h"
23#include "qgsmeshlayer.h"
26#include "qgsmeshlayerutils.h"
27#include "qgspolygon.h"
29#include "qgsrasterfilewriter.h"
30
31#include <QTextStream>
32
34
35
36static QgsFields createFields( const QList<QgsMeshDatasetGroupMetadata> &groupMetadataList, int vectorOption )
37{
38 QgsFields fields;
39 for ( const QgsMeshDatasetGroupMetadata &meta : groupMetadataList )
40 {
41 if ( meta.isVector() )
42 {
43 if ( vectorOption == 0 || vectorOption == 2 )
44 {
45 fields.append( QgsField( QStringLiteral( "%1_x" ).arg( meta.name() ), QMetaType::Type::Double ) );
46 fields.append( QgsField( QStringLiteral( "%1_y" ).arg( meta.name() ), QMetaType::Type::Double ) );
47 }
48
49 if ( vectorOption == 1 || vectorOption == 2 )
50 {
51 fields.append( QgsField( QStringLiteral( "%1_mag" ).arg( meta.name() ), QMetaType::Type::Double ) );
52 fields.append( QgsField( QStringLiteral( "%1_dir" ).arg( meta.name() ), QMetaType::Type::Double ) );
53 }
54 }
55 else
56 fields.append( QgsField( meta.name(), QMetaType::Type::Double ) );
57 }
58 return fields;
59}
60
61static QVector<double> vectorValue( const QgsMeshDatasetValue &value, int exportOption )
62{
63 QVector<double> ret( exportOption == 2 ? 4 : 2 );
64
65 if ( exportOption == 0 || exportOption == 2 )
66 {
67 ret[0] = value.x();
68 ret[1] = value.y();
69 }
70 if ( exportOption == 1 || exportOption == 2 )
71 {
72 double x = value.x();
73 double y = value.y();
74 double magnitude = sqrt( x * x + y * y );
75 double direction = ( asin( x / magnitude ) ) / M_PI * 180;
76 if ( y < 0 )
77 direction = 180 - direction;
78
79 if ( exportOption == 1 )
80 {
81 ret[0] = magnitude;
82 ret[1] = direction;
83 }
84 if ( exportOption == 2 )
85 {
86 ret[2] = magnitude;
87 ret[3] = direction;
88 }
89 }
90 return ret;
91}
92
93static void addAttributes( const QgsMeshDatasetValue &value, QgsAttributes &attributes, bool isVector, int vectorOption )
94{
95 if ( isVector )
96 {
97 QVector<double> vectorValues = vectorValue( value, vectorOption );
98 for ( double v : vectorValues )
99 {
100 if ( v == std::numeric_limits<double>::quiet_NaN() )
101 attributes.append( QVariant() );
102 else
103 attributes.append( v );
104 }
105 }
106 else
107 {
108 if ( value.scalar() == std::numeric_limits<double>::quiet_NaN() )
109 attributes.append( QVariant() );
110 else
111 attributes.append( value.scalar() );
112 }
113}
114
115static QgsMeshDatasetValue extractDatasetValue(
116 const QgsPointXY &point,
117 int nativeFaceIndex,
118 int triangularFaceIndex,
119 const QgsTriangularMesh &triangularMesh,
120 const QgsMeshDataBlock &activeFaces,
121 const QgsMeshDataBlock &datasetValues,
122 const QgsMeshDatasetGroupMetadata &metadata
123)
124{
125 bool faceActive = activeFaces.active( nativeFaceIndex );
127 if ( faceActive )
128 {
129 switch ( metadata.dataType() )
130 {
132 //not supported
133 break;
136 {
137 value = datasetValues.value( nativeFaceIndex );
138 }
139 break;
140
142 {
143 const QgsMeshFace &face = triangularMesh.triangles()[triangularFaceIndex];
144 const int v1 = face[0], v2 = face[1], v3 = face[2];
145 const QgsPoint p1 = triangularMesh.vertices()[v1], p2 = triangularMesh.vertices()[v2], p3 = triangularMesh.vertices()[v3];
146 const QgsMeshDatasetValue val1 = datasetValues.value( v1 );
147 const QgsMeshDatasetValue val2 = datasetValues.value( v2 );
148 const QgsMeshDatasetValue val3 = datasetValues.value( v3 );
149 const double x = QgsMeshLayerUtils::interpolateFromVerticesData( p1, p2, p3, val1.x(), val2.x(), val3.x(), point );
150 double y = std::numeric_limits<double>::quiet_NaN();
151 bool isVector = metadata.isVector();
152 if ( isVector )
153 y = QgsMeshLayerUtils::interpolateFromVerticesData( p1, p2, p3, val1.y(), val2.y(), val3.y(), point );
154
155 value = QgsMeshDatasetValue( x, y );
156 }
157 break;
158 }
159 }
160
161 return value;
162}
163
164QString QgsExportMeshOnElement::group() const
165{
166 return QObject::tr( "Mesh" );
167}
168
169QString QgsExportMeshOnElement::groupId() const
170{
171 return QStringLiteral( "mesh" );
172}
173
174QString QgsExportMeshVerticesAlgorithm::shortHelpString() const
175{
176 return QObject::tr( "This algorithm exports a mesh layer's vertices to a point vector layer, with the dataset values on vertices as attribute values." );
177}
178
179QString QgsExportMeshVerticesAlgorithm::shortDescription() const
180{
181 return QObject::tr( "Exports mesh vertices to a point vector layer." );
182}
183
184QString QgsExportMeshVerticesAlgorithm::name() const
185{
186 return QStringLiteral( "exportmeshvertices" );
187}
188
189QString QgsExportMeshVerticesAlgorithm::displayName() const
190{
191 return QObject::tr( "Export mesh vertices" );
192}
193
194QgsProcessingAlgorithm *QgsExportMeshVerticesAlgorithm::createInstance() const
195{
196 return new QgsExportMeshVerticesAlgorithm();
197}
198
199QgsGeometry QgsExportMeshVerticesAlgorithm::meshElement( int index ) const
200{
201 return QgsGeometry( new QgsPoint( mNativeMesh.vertex( index ) ) );
202}
203
204void QgsExportMeshOnElement::initAlgorithm( const QVariantMap &configuration )
205{
206 Q_UNUSED( configuration );
207
208 addParameter( new QgsProcessingParameterMeshLayer( QStringLiteral( "INPUT" ), QObject::tr( "Input mesh layer" ) ) );
209
210
212 QStringLiteral( "DATASET_GROUPS" ),
213 QObject::tr( "Dataset groups" ),
214 QStringLiteral( "INPUT" ),
215 supportedDataType(), true
216 ) );
217
219 QStringLiteral( "DATASET_TIME" ),
220 QObject::tr( "Dataset time" ),
221 QStringLiteral( "INPUT" ),
222 QStringLiteral( "DATASET_GROUPS" )
223 ) );
224
225 addParameter( new QgsProcessingParameterCrs( QStringLiteral( "CRS_OUTPUT" ), QObject::tr( "Output coordinate system" ), QVariant(), true ) );
226
227 QStringList exportVectorOptions;
228 exportVectorOptions << QObject::tr( "Cartesian (x,y)" )
229 << QObject::tr( "Polar (magnitude,degree)" )
230 << QObject::tr( "Cartesian and Polar" );
231 addParameter( new QgsProcessingParameterEnum( QStringLiteral( "VECTOR_OPTION" ), QObject::tr( "Export vector option" ), exportVectorOptions, false, 0 ) );
232 addParameter( new QgsProcessingParameterFeatureSink( QStringLiteral( "OUTPUT" ), QObject::tr( "Output vector layer" ), sinkType() ) );
233}
234
235static QgsInterval datasetRelativetime( const QVariant parameterTimeVariant, QgsMeshLayer *meshLayer, const QgsProcessingContext &context )
236{
237 QgsInterval relativeTime( 0 );
238 QDateTime layerReferenceTime = static_cast<QgsMeshLayerTemporalProperties *>( meshLayer->temporalProperties() )->referenceTime();
239 QString timeType = QgsProcessingParameterMeshDatasetTime::valueAsTimeType( parameterTimeVariant );
240
241 if ( timeType == QLatin1String( "dataset-time-step" ) )
242 {
244 relativeTime = meshLayer->datasetRelativeTime( datasetIndex );
245 }
246 else if ( timeType == QLatin1String( "defined-date-time" ) )
247 {
248 QDateTime dateTime = QgsProcessingParameterMeshDatasetTime::timeValueAsDefinedDateTime( parameterTimeVariant );
249 if ( dateTime.isValid() )
250 relativeTime = QgsInterval( layerReferenceTime.secsTo( dateTime ) );
251 }
252 else if ( timeType == QLatin1String( "current-context-time" ) )
253 {
254 QDateTime dateTime = context.currentTimeRange().begin();
255 if ( dateTime.isValid() )
256 relativeTime = QgsInterval( layerReferenceTime.secsTo( dateTime ) );
257 }
258
259 return relativeTime;
260}
261
262
263bool QgsExportMeshOnElement::prepareAlgorithm( const QVariantMap &parameters, QgsProcessingContext &context, QgsProcessingFeedback *feedback )
264{
265 QgsMeshLayer *meshLayer = parameterAsMeshLayer( parameters, QStringLiteral( "INPUT" ), context );
266
267 if ( !meshLayer || !meshLayer->isValid() )
268 return false;
269
270 if ( meshLayer->isEditable() )
271 throw QgsProcessingException( QObject::tr( "Input mesh layer in edit mode is not supported" ) );
272
273 QgsCoordinateReferenceSystem outputCrs = parameterAsCrs( parameters, QStringLiteral( "CRS_OUTPUT" ), context );
274 if ( !outputCrs.isValid() )
275 outputCrs = meshLayer->crs();
276 mTransform = QgsCoordinateTransform( meshLayer->crs(), outputCrs, context.transformContext() );
277 if ( !meshLayer->nativeMesh() )
278 meshLayer->updateTriangularMesh( mTransform ); //necessary to load the native mesh
279
280 mNativeMesh = *meshLayer->nativeMesh();
281
282 QList<int> datasetGroups = QgsProcessingParameterMeshDatasetGroups::valueAsDatasetGroup( parameters.value( QStringLiteral( "DATASET_GROUPS" ) ) );
283
284 if ( feedback )
285 {
286 feedback->setProgressText( QObject::tr( "Preparing data" ) );
287 }
288
289 // Extract the date time used to export dataset values under a relative time
290 QVariant parameterTimeVariant = parameters.value( QStringLiteral( "DATASET_TIME" ) );
291 QgsInterval relativeTime = datasetRelativetime( parameterTimeVariant, meshLayer, context );
292
293 switch ( meshElementType() )
294 {
295 case QgsMesh::Face:
296 mElementCount = mNativeMesh.faceCount();
297 break;
298 case QgsMesh::Vertex:
299 mElementCount = mNativeMesh.vertexCount();
300 break;
301 case QgsMesh::Edge:
302 mElementCount = mNativeMesh.edgeCount();
303 break;
304 }
305
306 for ( int i = 0; i < datasetGroups.count(); ++i )
307 {
308 int groupIndex = datasetGroups.at( i );
309 QgsMeshDatasetIndex datasetIndex = meshLayer->datasetIndexAtRelativeTime( relativeTime, groupIndex );
310
311 DataGroup dataGroup;
312 dataGroup.metadata = meshLayer->datasetGroupMetadata( datasetIndex );
313 if ( supportedDataType().contains( dataGroup.metadata.dataType() ) )
314 {
315 dataGroup.datasetValues = meshLayer->datasetValues( datasetIndex, 0, mElementCount );
316 mDataPerGroup.append( dataGroup );
317 }
318 if ( feedback )
319 feedback->setProgress( 100 * i / datasetGroups.count() );
320 }
321
322 mExportVectorOption = parameterAsInt( parameters, QStringLiteral( "VECTOR_OPTION" ), context );
323
324 return true;
325}
326
327QVariantMap QgsExportMeshOnElement::processAlgorithm( const QVariantMap &parameters, QgsProcessingContext &context, QgsProcessingFeedback *feedback )
328{
329 if ( feedback )
330 {
331 if ( feedback->isCanceled() )
332 return QVariantMap();
333 feedback->setProgress( 0 );
334 feedback->setProgressText( QObject::tr( "Creating output vector layer" ) );
335 }
336
337 QList<QgsMeshDatasetGroupMetadata> metaList;
338 metaList.reserve( mDataPerGroup.size() );
339 for ( const DataGroup &dataGroup : std::as_const( mDataPerGroup ) )
340 metaList.append( dataGroup.metadata );
341 QgsFields fields = createFields( metaList, mExportVectorOption );
342
343 QgsCoordinateReferenceSystem outputCrs = parameterAsCrs( parameters, QStringLiteral( "CRS_OUTPUT" ), context );
344 QString identifier;
345 std::unique_ptr<QgsFeatureSink> sink( parameterAsSink( parameters, QStringLiteral( "OUTPUT" ), context, identifier, fields, sinkGeometryType(), outputCrs ) );
346 if ( !sink )
347 return QVariantMap();
348
349 if ( feedback )
350 {
351 if ( feedback->isCanceled() )
352 return QVariantMap();
353 feedback->setProgress( 0 );
354 feedback->setProgressText( QObject::tr( "Creating points for each vertices" ) );
355 }
356
357 for ( int i = 0; i < mElementCount; ++i )
358 {
359 QgsAttributes attributes;
360 for ( const DataGroup &dataGroup : std::as_const( mDataPerGroup ) )
361 {
362 const QgsMeshDatasetValue &value = dataGroup.datasetValues.value( i );
363 addAttributes( value, attributes, dataGroup.metadata.isVector(), mExportVectorOption );
364 }
365
366 QgsFeature feat;
367 QgsGeometry geom = meshElement( i );
368 try
369 {
370 geom.transform( mTransform );
371 }
372 catch ( QgsCsException & )
373 {
374 geom = meshElement( i );
375 if ( feedback )
376 feedback->reportError( QObject::tr( "Could not transform point to destination CRS" ) );
377 }
378 feat.setGeometry( geom );
379 feat.setAttributes( attributes );
380
381 if ( !sink->addFeature( feat, QgsFeatureSink::FastInsert ) )
382 throw QgsProcessingException( writeFeatureError( sink.get(), parameters, QStringLiteral( "OUTPUT" ) ) );
383
384 if ( feedback )
385 {
386 if ( feedback->isCanceled() )
387 return QVariantMap();
388 feedback->setProgress( 100 * i / mElementCount );
389 }
390 }
391
392 sink->finalize();
393
394 QVariantMap ret;
395 ret[QStringLiteral( "OUTPUT" )] = identifier;
396
397 return ret;
398}
399
400QString QgsExportMeshFacesAlgorithm::shortHelpString() const
401{
402 return QObject::tr( "This algorithm exports a mesh layer's faces to a polygon vector layer, with the dataset values on faces as attribute values." );
403}
404
405QString QgsExportMeshFacesAlgorithm::shortDescription() const
406{
407 return QObject::tr( "Exports mesh faces to a polygon vector layer." );
408}
409
410QString QgsExportMeshFacesAlgorithm::name() const
411{
412 return QStringLiteral( "exportmeshfaces" );
413}
414
415QString QgsExportMeshFacesAlgorithm::displayName() const
416{
417 return QObject::tr( "Export mesh faces" );
418}
419
420QgsProcessingAlgorithm *QgsExportMeshFacesAlgorithm::createInstance() const
421{
422 return new QgsExportMeshFacesAlgorithm();
423}
424
425QgsGeometry QgsExportMeshFacesAlgorithm::meshElement( int index ) const
426{
427 const QgsMeshFace &face = mNativeMesh.face( index );
428 QVector<QgsPoint> vertices( face.size() );
429 for ( int i = 0; i < face.size(); ++i )
430 vertices[i] = mNativeMesh.vertex( face.at( i ) );
431 auto polygon = std::make_unique<QgsPolygon>();
432 polygon->setExteriorRing( new QgsLineString( vertices ) );
433 return QgsGeometry( polygon.release() );
434}
435
436QString QgsExportMeshEdgesAlgorithm::shortHelpString() const
437{
438 return QObject::tr( "This algorithm exports a mesh layer's edges to a line vector layer, with the dataset values on edges as attribute values." );
439}
440
441QString QgsExportMeshEdgesAlgorithm::shortDescription() const
442{
443 return QObject::tr( "Exports mesh edges to a line vector layer." );
444}
445
446QString QgsExportMeshEdgesAlgorithm::name() const
447{
448 return QStringLiteral( "exportmeshedges" );
449}
450
451QString QgsExportMeshEdgesAlgorithm::displayName() const
452{
453 return QObject::tr( "Export mesh edges" );
454}
455
456QgsProcessingAlgorithm *QgsExportMeshEdgesAlgorithm::createInstance() const
457{
458 return new QgsExportMeshEdgesAlgorithm();
459}
460
461QgsGeometry QgsExportMeshEdgesAlgorithm::meshElement( int index ) const
462{
463 const QgsMeshEdge &edge = mNativeMesh.edge( index );
464 QVector<QgsPoint> vertices( 2 );
465 vertices[0] = mNativeMesh.vertex( edge.first );
466 vertices[1] = mNativeMesh.vertex( edge.second );
467 return QgsGeometry( new QgsLineString( vertices ) );
468}
469
470
471QString QgsExportMeshOnGridAlgorithm::name() const { return QStringLiteral( "exportmeshongrid" ); }
472
473QString QgsExportMeshOnGridAlgorithm::displayName() const { return QObject::tr( "Export mesh on grid" ); }
474
475QString QgsExportMeshOnGridAlgorithm::group() const { return QObject::tr( "Mesh" ); }
476
477QString QgsExportMeshOnGridAlgorithm::groupId() const { return QStringLiteral( "mesh" ); }
478
479QString QgsExportMeshOnGridAlgorithm::shortHelpString() const
480{
481 return QObject::tr( "This algorithm exports a mesh layer's dataset values to a gridded point vector layer, with the dataset values on each point as attribute values.\n"
482 "For data on volume (3D stacked dataset values), the exported dataset values are averaged on faces using the method defined in the mesh layer properties (default is Multi level averaging method).\n"
483 "1D meshes are not supported." );
484}
485
486QString QgsExportMeshOnGridAlgorithm::shortDescription() const
487{
488 return QObject::tr( "Exports mesh dataset values to a gridded point vector layer." );
489}
490
491QgsProcessingAlgorithm *QgsExportMeshOnGridAlgorithm::createInstance() const
492{
493 return new QgsExportMeshOnGridAlgorithm();
494}
495
496void QgsExportMeshOnGridAlgorithm::initAlgorithm( const QVariantMap &configuration )
497{
498 Q_UNUSED( configuration );
499
500 addParameter( new QgsProcessingParameterMeshLayer( QStringLiteral( "INPUT" ), QObject::tr( "Input mesh layer" ) ) );
501
503 QStringLiteral( "DATASET_GROUPS" ),
504 QObject::tr( "Dataset groups" ),
505 QStringLiteral( "INPUT" ),
506 supportedDataType()
507 ) );
508
510 QStringLiteral( "DATASET_TIME" ),
511 QObject::tr( "Dataset time" ),
512 QStringLiteral( "INPUT" ),
513 QStringLiteral( "DATASET_GROUPS" )
514 ) );
515
516 addParameter( new QgsProcessingParameterExtent( QStringLiteral( "EXTENT" ), QObject::tr( "Extent" ), QVariant(), true ) );
517
518 addParameter( new QgsProcessingParameterDistance( QStringLiteral( "GRID_SPACING" ), QObject::tr( "Grid spacing" ), 10, QStringLiteral( "INPUT" ), false ) );
519
520 addParameter( new QgsProcessingParameterCrs( QStringLiteral( "CRS_OUTPUT" ), QObject::tr( "Output coordinate system" ), QVariant(), true ) );
521
522 QStringList exportVectorOptions;
523 exportVectorOptions << QObject::tr( "Cartesian (x,y)" )
524 << QObject::tr( "Polar (magnitude,degree)" )
525 << QObject::tr( "Cartesian and Polar" );
526 addParameter( new QgsProcessingParameterEnum( QStringLiteral( "VECTOR_OPTION" ), QObject::tr( "Export vector option" ), exportVectorOptions, false, 0 ) );
527 addParameter( new QgsProcessingParameterFeatureSink( QStringLiteral( "OUTPUT" ), QObject::tr( "Output vector layer" ), Qgis::ProcessingSourceType::VectorPoint ) );
528}
529
530static void extractDatasetValues( const QList<int> &datasetGroups, QgsMeshLayer *meshLayer, const QgsMesh &nativeMesh, const QgsInterval &relativeTime, const QSet<int> supportedDataType, QList<DataGroup> &datasetPerGroup, QgsProcessingFeedback *feedback )
531{
532 for ( int i = 0; i < datasetGroups.count(); ++i )
533 {
534 int groupIndex = datasetGroups.at( i );
535 QgsMeshDatasetIndex datasetIndex = meshLayer->datasetIndexAtRelativeTime( relativeTime, groupIndex );
536
537 DataGroup dataGroup;
538 dataGroup.metadata = meshLayer->datasetGroupMetadata( datasetIndex );
539 if ( supportedDataType.contains( dataGroup.metadata.dataType() ) )
540 {
541 int valueCount = dataGroup.metadata.dataType() == QgsMeshDatasetGroupMetadata::DataOnVertices ? nativeMesh.vertices.count() : nativeMesh.faceCount();
542 dataGroup.datasetValues = meshLayer->datasetValues( datasetIndex, 0, valueCount );
543 dataGroup.activeFaces = meshLayer->areFacesActive( datasetIndex, 0, nativeMesh.faceCount() );
544 if ( dataGroup.metadata.dataType() == QgsMeshDatasetGroupMetadata::DataOnVolumes )
545 {
546 dataGroup.dataset3dStakedValue = meshLayer->dataset3dValues( datasetIndex, 0, valueCount );
547 }
548 datasetPerGroup.append( dataGroup );
549 }
550 if ( feedback )
551 feedback->setProgress( 100 * i / datasetGroups.count() );
552 }
553}
554
555bool QgsExportMeshOnGridAlgorithm::prepareAlgorithm( const QVariantMap &parameters, QgsProcessingContext &context, QgsProcessingFeedback *feedback )
556{
557 QgsMeshLayer *meshLayer = parameterAsMeshLayer( parameters, QStringLiteral( "INPUT" ), context );
558
559 if ( !meshLayer || !meshLayer->isValid() )
560 return false;
561
562 QgsCoordinateReferenceSystem outputCrs = parameterAsCrs( parameters, QStringLiteral( "CRS_OUTPUT" ), context );
563 if ( !outputCrs.isValid() )
564 outputCrs = meshLayer->crs();
565 mTransform = QgsCoordinateTransform( meshLayer->crs(), outputCrs, context.transformContext() );
566 if ( !meshLayer->nativeMesh() )
567 meshLayer->updateTriangularMesh( mTransform ); //necessary to load the native mesh
568
569 const QgsMesh &nativeMesh = *meshLayer->nativeMesh();
570
571 QList<int> datasetGroups = QgsProcessingParameterMeshDatasetGroups::valueAsDatasetGroup( parameters.value( QStringLiteral( "DATASET_GROUPS" ) ) );
572
573 if ( feedback )
574 {
575 feedback->setProgressText( QObject::tr( "Preparing data" ) );
576 }
577
578 // Extract the date time used to export dataset values under a relative time
579 QVariant parameterTimeVariant = parameters.value( QStringLiteral( "DATASET_TIME" ) );
580 QgsInterval relativeTime = datasetRelativetime( parameterTimeVariant, meshLayer, context );
581
582 extractDatasetValues( datasetGroups, meshLayer, nativeMesh, relativeTime, supportedDataType(), mDataPerGroup, feedback );
583 mTriangularMesh.update( meshLayer->nativeMesh(), mTransform );
584
585 mExportVectorOption = parameterAsInt( parameters, QStringLiteral( "VECTOR_OPTION" ), context );
586
587 return true;
588}
589
590QVariantMap QgsExportMeshOnGridAlgorithm::processAlgorithm( const QVariantMap &parameters, QgsProcessingContext &context, QgsProcessingFeedback *feedback )
591{
592 if ( feedback )
593 {
594 if ( feedback->isCanceled() )
595 return QVariantMap();
596 feedback->setProgress( 0 );
597 feedback->setProgressText( QObject::tr( "Creating output vector layer" ) );
598 }
599
600 //First, if present, average 3D staked dataset value to 2D face value
601 const QgsMesh3DAveragingMethod *avgMethod = mLayerRendererSettings.averagingMethod();
602 for ( DataGroup &dataGroup : mDataPerGroup )
603 {
604 if ( dataGroup.dataset3dStakedValue.isValid() )
605 dataGroup.datasetValues = avgMethod->calculate( dataGroup.dataset3dStakedValue );
606 }
607
608 QList<QgsMeshDatasetGroupMetadata> metaList;
609 metaList.reserve( mDataPerGroup.size() );
610 for ( const DataGroup &dataGroup : std::as_const( mDataPerGroup ) )
611 metaList.append( dataGroup.metadata );
612 QgsFields fields = createFields( metaList, mExportVectorOption );
613
614 //create sink
615 QgsCoordinateReferenceSystem outputCrs = parameterAsCrs( parameters, QStringLiteral( "CRS_OUTPUT" ), context );
616 QString identifier;
617 std::unique_ptr<QgsFeatureSink> sink( parameterAsSink( parameters, QStringLiteral( "OUTPUT" ), context, identifier, fields, Qgis::WkbType::Point, outputCrs ) );
618 if ( !sink )
619 return QVariantMap();
620
621 if ( feedback )
622 {
623 if ( feedback->isCanceled() )
624 return QVariantMap();
625 feedback->setProgress( 0 );
626 feedback->setProgressText( QObject::tr( "Creating gridded points" ) );
627 }
628
629 // grid definition
630 const double gridSpacing = parameterAsDouble( parameters, QStringLiteral( "GRID_SPACING" ), context );
631 if ( qgsDoubleNear( gridSpacing, 0 ) )
632 {
633 throw QgsProcessingException( QObject::tr( "Grid spacing cannot be 0" ) );
634 }
635
636 QgsRectangle extent = parameterAsExtent( parameters, QStringLiteral( "EXTENT" ), context );
637 if ( extent.isEmpty() )
638 extent = mTriangularMesh.extent();
639 int pointXCount = int( extent.width() / gridSpacing ) + 1;
640 int pointYCount = int( extent.height() / gridSpacing ) + 1;
641
642 for ( int ix = 0; ix < pointXCount; ++ix )
643 {
644 for ( int iy = 0; iy < pointYCount; ++iy )
645 {
646 QgsPoint point( extent.xMinimum() + ix * gridSpacing, extent.yMinimum() + iy * gridSpacing );
647 int triangularFaceIndex = mTriangularMesh.faceIndexForPoint_v2( point );
648 if ( triangularFaceIndex >= 0 )
649 {
650 //extract dataset values for the point
651 QgsAttributes attributes;
652 int nativeFaceIndex = mTriangularMesh.trianglesToNativeFaces().at( triangularFaceIndex );
653 for ( int i = 0; i < mDataPerGroup.count(); ++i )
654 {
655 const DataGroup &dataGroup = mDataPerGroup.at( i );
656 bool faceActive = dataGroup.activeFaces.active( nativeFaceIndex );
657 if ( !faceActive )
658 continue;
659 QgsMeshDatasetValue value = extractDatasetValue(
660 point,
661 nativeFaceIndex,
662 triangularFaceIndex,
663 mTriangularMesh,
664 dataGroup.activeFaces,
665 dataGroup.datasetValues,
666 dataGroup.metadata
667 );
668
669 if ( dataGroup.metadata.isVector() )
670 {
671 QVector<double> vector = vectorValue( dataGroup.datasetValues.value( i ), mExportVectorOption );
672 for ( double v : vector )
673 {
674 attributes.append( v );
675 }
676 }
677 else
678 attributes.append( value.scalar() );
679 }
680 QgsFeature feat;
681 QgsGeometry geom( point.clone() );
682 try
683 {
684 geom.transform( mTransform );
685 }
686 catch ( QgsCsException & )
687 {
688 geom = QgsGeometry( point.clone() );
689 feedback->reportError( QObject::tr( "Could not transform point to destination CRS" ) );
690 }
691 feat.setGeometry( geom );
692 feat.setAttributes( attributes );
693
694 if ( !sink->addFeature( feat, QgsFeatureSink::FastInsert ) )
695 {
696 throw QgsProcessingException( writeFeatureError( sink.get(), parameters, QString() ) );
697 }
698 }
699 }
700 }
701
702 sink->finalize();
703
704 QVariantMap ret;
705 ret[QStringLiteral( "OUTPUT" )] = identifier;
706
707 return ret;
708}
709
710QSet<int> QgsExportMeshOnGridAlgorithm::supportedDataType()
711{
712 return QSet<int>(
716 );
717}
718
719QString QgsMeshRasterizeAlgorithm::name() const
720{
721 return QStringLiteral( "meshrasterize" );
722}
723
724QString QgsMeshRasterizeAlgorithm::displayName() const
725{
726 return QObject::tr( "Rasterize mesh dataset" );
727}
728
729QString QgsMeshRasterizeAlgorithm::group() const
730{
731 return QObject::tr( "Mesh" );
732}
733
734QString QgsMeshRasterizeAlgorithm::groupId() const
735{
736 return QStringLiteral( "mesh" );
737}
738
739QString QgsMeshRasterizeAlgorithm::shortHelpString() const
740{
741 return QObject::tr( "This algorithm creates a raster layer from a mesh dataset.\n"
742 "For data on volume (3D stacked dataset values), the exported dataset values are averaged on faces using the method defined in the mesh layer properties (default is Multi level averaging method).\n"
743 "1D meshes are not supported." );
744}
745
746QString QgsMeshRasterizeAlgorithm::shortDescription() const
747{
748 return QObject::tr( "Creates a raster layer from a mesh dataset." );
749}
750
751QgsProcessingAlgorithm *QgsMeshRasterizeAlgorithm::createInstance() const
752{
753 return new QgsMeshRasterizeAlgorithm();
754}
755
756void QgsMeshRasterizeAlgorithm::initAlgorithm( const QVariantMap &configuration )
757{
758 Q_UNUSED( configuration );
759
760 addParameter( new QgsProcessingParameterMeshLayer( QStringLiteral( "INPUT" ), QObject::tr( "Input mesh layer" ) ) );
761
763 QStringLiteral( "DATASET_GROUPS" ),
764 QObject::tr( "Dataset groups" ),
765 QStringLiteral( "INPUT" ),
766 supportedDataType(),
767 true
768 ) );
769
771 QStringLiteral( "DATASET_TIME" ),
772 QObject::tr( "Dataset time" ),
773 QStringLiteral( "INPUT" ),
774 QStringLiteral( "DATASET_GROUPS" )
775 ) );
776
777 addParameter( new QgsProcessingParameterExtent( QStringLiteral( "EXTENT" ), QObject::tr( "Extent" ), QVariant(), true ) );
778 addParameter( new QgsProcessingParameterDistance( QStringLiteral( "PIXEL_SIZE" ), QObject::tr( "Pixel size" ), 1, QStringLiteral( "INPUT" ), false ) );
779 addParameter( new QgsProcessingParameterCrs( QStringLiteral( "CRS_OUTPUT" ), QObject::tr( "Output coordinate system" ), QVariant(), true ) );
780
781 // backwards compatibility parameter
782 // TODO QGIS 4: remove parameter and related logic
783 auto createOptsParam = std::make_unique<QgsProcessingParameterString>( QStringLiteral( "CREATE_OPTIONS" ), QObject::tr( "Creation options" ), QVariant(), false, true );
784 createOptsParam->setMetadata( QVariantMap( { { QStringLiteral( "widget_wrapper" ), QVariantMap( { { QStringLiteral( "widget_type" ), QStringLiteral( "rasteroptions" ) } } ) } } ) );
785 createOptsParam->setFlags( createOptsParam->flags() | Qgis::ProcessingParameterFlag::Hidden );
786 addParameter( createOptsParam.release() );
787
788 auto creationOptsParam = std::make_unique<QgsProcessingParameterString>( QStringLiteral( "CREATION_OPTIONS" ), QObject::tr( "Creation options" ), QVariant(), false, true );
789 creationOptsParam->setMetadata( QVariantMap( { { QStringLiteral( "widget_wrapper" ), QVariantMap( { { QStringLiteral( "widget_type" ), QStringLiteral( "rasteroptions" ) } } ) } } ) );
790 creationOptsParam->setFlags( creationOptsParam->flags() | Qgis::ProcessingParameterFlag::Advanced );
791 addParameter( creationOptsParam.release() );
792
793 addParameter( new QgsProcessingParameterRasterDestination( QStringLiteral( "OUTPUT" ), QObject::tr( "Output raster layer" ) ) );
794}
795
796bool QgsMeshRasterizeAlgorithm::prepareAlgorithm( const QVariantMap &parameters, QgsProcessingContext &context, QgsProcessingFeedback *feedback )
797{
798 QgsMeshLayer *meshLayer = parameterAsMeshLayer( parameters, QStringLiteral( "INPUT" ), context );
799
800 if ( !meshLayer || !meshLayer->isValid() )
801 return false;
802
803 QgsCoordinateReferenceSystem outputCrs = parameterAsCrs( parameters, QStringLiteral( "CRS_OUTPUT" ), context );
804 if ( !outputCrs.isValid() )
805 outputCrs = meshLayer->crs();
806 mTransform = QgsCoordinateTransform( meshLayer->crs(), outputCrs, context.transformContext() );
807 if ( !meshLayer->nativeMesh() )
808 meshLayer->updateTriangularMesh( mTransform ); //necessary to load the native mesh
809
810 mTriangularMesh.update( meshLayer->nativeMesh(), mTransform );
811
812 QList<int> datasetGroups = QgsProcessingParameterMeshDatasetGroups::valueAsDatasetGroup( parameters.value( QStringLiteral( "DATASET_GROUPS" ) ) );
813
814 if ( feedback )
815 {
816 feedback->setProgressText( QObject::tr( "Preparing data" ) );
817 }
818
819 // Extract the date time used to export dataset values under a relative time
820 QVariant parameterTimeVariant = parameters.value( QStringLiteral( "DATASET_TIME" ) );
821 QgsInterval relativeTime = datasetRelativetime( parameterTimeVariant, meshLayer, context );
822
823 extractDatasetValues( datasetGroups, meshLayer, *meshLayer->nativeMesh(), relativeTime, supportedDataType(), mDataPerGroup, feedback );
824
825 mLayerRendererSettings = meshLayer->rendererSettings();
826
827 return true;
828}
829
830QVariantMap QgsMeshRasterizeAlgorithm::processAlgorithm( const QVariantMap &parameters, QgsProcessingContext &context, QgsProcessingFeedback *feedback )
831{
832 if ( feedback )
833 {
834 if ( feedback->isCanceled() )
835 return QVariantMap();
836 feedback->setProgress( 0 );
837 feedback->setProgressText( QObject::tr( "Creating raster layer" ) );
838 }
839
840 //First, if present, average 3D staked dataset value to 2D face value
841 const QgsMesh3DAveragingMethod *avgMethod = mLayerRendererSettings.averagingMethod();
842 for ( DataGroup &dataGroup : mDataPerGroup )
843 {
844 if ( dataGroup.dataset3dStakedValue.isValid() )
845 dataGroup.datasetValues = avgMethod->calculate( dataGroup.dataset3dStakedValue );
846 }
847
848 // create raster
849 const double pixelSize = parameterAsDouble( parameters, QStringLiteral( "PIXEL_SIZE" ), context );
850 if ( qgsDoubleNear( pixelSize, 0 ) )
851 {
852 throw QgsProcessingException( QObject::tr( "Pixel size cannot be 0" ) );
853 }
854
855 QgsRectangle extent = parameterAsExtent( parameters, QStringLiteral( "EXTENT" ), context );
856 if ( extent.isEmpty() )
857 extent = mTriangularMesh.extent();
858
859 int width = extent.width() / pixelSize;
860 int height = extent.height() / pixelSize;
861
862 QString creationOptions = parameterAsString( parameters, QStringLiteral( "CREATION_OPTIONS" ), context ).trimmed();
863 // handle backwards compatibility parameter CREATE_OPTIONS
864 const QString optionsString = parameterAsString( parameters, QStringLiteral( "CREATE_OPTIONS" ), context );
865 if ( !optionsString.isEmpty() )
866 creationOptions = optionsString;
867
868 const QString fileName = parameterAsOutputLayer( parameters, QStringLiteral( "OUTPUT" ), context );
869 const QFileInfo fileInfo( fileName );
870 const QString outputFormat = QgsRasterFileWriter::driverForExtension( fileInfo.suffix() );
871 QgsRasterFileWriter rasterFileWriter( fileName );
872 rasterFileWriter.setOutputProviderKey( QStringLiteral( "gdal" ) );
873 if ( !creationOptions.isEmpty() )
874 {
875 rasterFileWriter.setCreationOptions( creationOptions.split( '|' ) );
876 }
877 rasterFileWriter.setOutputFormat( outputFormat );
878
879 std::unique_ptr<QgsRasterDataProvider> rasterDataProvider(
880 rasterFileWriter.createMultiBandRaster( Qgis::DataType::Float64, width, height, extent, mTransform.destinationCrs(), mDataPerGroup.count() )
881 );
882 rasterDataProvider->setEditable( true );
883
884 for ( int i = 0; i < mDataPerGroup.count(); ++i )
885 {
886 const DataGroup &dataGroup = mDataPerGroup.at( i );
887 QgsRasterBlockFeedback rasterBlockFeedBack;
888 if ( feedback )
889 QObject::connect( &rasterBlockFeedBack, &QgsFeedback::canceled, feedback, &QgsFeedback::cancel );
890
891 if ( dataGroup.datasetValues.isValid() )
892 {
893 std::unique_ptr<QgsRasterBlock> block( QgsMeshUtils::exportRasterBlock(
894 mTriangularMesh,
895 dataGroup.datasetValues,
896 dataGroup.activeFaces,
897 dataGroup.metadata.dataType(),
898 mTransform,
899 pixelSize,
900 extent,
901 &rasterBlockFeedBack
902 ) );
903
904 if ( !rasterDataProvider->writeBlock( block.get(), i + 1 ) )
905 {
906 throw QgsProcessingException( QObject::tr( "Could not write raster block: %1" ).arg( rasterDataProvider->error().summary() ) );
907 }
908 rasterDataProvider->setNoDataValue( i + 1, block->noDataValue() );
909 }
910 else
911 rasterDataProvider->setNoDataValue( i + 1, std::numeric_limits<double>::quiet_NaN() );
912
913 if ( feedback )
914 {
915 if ( feedback->isCanceled() )
916 return QVariantMap();
917 feedback->setProgress( 100 * i / mDataPerGroup.count() );
918 }
919 }
920
921 rasterDataProvider->setEditable( false );
922
923 if ( feedback )
924 feedback->setProgress( 100 );
925
926 QVariantMap ret;
927 ret[QStringLiteral( "OUTPUT" )] = fileName;
928
929 return ret;
930}
931
932QSet<int> QgsMeshRasterizeAlgorithm::supportedDataType()
933{
934 return QSet<int>(
938 );
939}
940
941QString QgsMeshContoursAlgorithm::name() const
942{
943 return QStringLiteral( "meshcontours" );
944}
945
946QString QgsMeshContoursAlgorithm::displayName() const
947{
948 return QObject::tr( "Export contours" );
949}
950
951QString QgsMeshContoursAlgorithm::group() const
952{
953 return QObject::tr( "Mesh" );
954}
955
956QString QgsMeshContoursAlgorithm::groupId() const
957{
958 return QStringLiteral( "mesh" );
959}
960
961QString QgsMeshContoursAlgorithm::shortHelpString() const
962{
963 return QObject::tr( "This algorithm creates contours as a vector layer from a mesh scalar dataset." );
964}
965
966QString QgsMeshContoursAlgorithm::shortDescription() const
967{
968 return QObject::tr( "Creates contours as vector layer from mesh scalar dataset." );
969}
970
971QgsProcessingAlgorithm *QgsMeshContoursAlgorithm::createInstance() const
972{
973 return new QgsMeshContoursAlgorithm();
974}
975
976void QgsMeshContoursAlgorithm::initAlgorithm( const QVariantMap &configuration )
977{
978 Q_UNUSED( configuration );
979
980 addParameter( new QgsProcessingParameterMeshLayer( QStringLiteral( "INPUT" ), QObject::tr( "Input mesh layer" ) ) );
981
983 QStringLiteral( "DATASET_GROUPS" ),
984 QObject::tr( "Dataset groups" ),
985 QStringLiteral( "INPUT" ),
986 supportedDataType()
987 ) );
988
990 QStringLiteral( "DATASET_TIME" ),
991 QObject::tr( "Dataset time" ),
992 QStringLiteral( "INPUT" ),
993 QStringLiteral( "DATASET_GROUPS" )
994 ) );
995
996 addParameter( new QgsProcessingParameterNumber(
997 QStringLiteral( "INCREMENT" ), QObject::tr( "Increment between contour levels" ), Qgis::ProcessingNumberParameterType::Double, QVariant(), true
998 ) );
999
1000 addParameter( new QgsProcessingParameterNumber(
1001 QStringLiteral( "MINIMUM" ), QObject::tr( "Minimum contour level" ), Qgis::ProcessingNumberParameterType::Double, QVariant(), true
1002 ) );
1003 addParameter( new QgsProcessingParameterNumber(
1004 QStringLiteral( "MAXIMUM" ), QObject::tr( "Maximum contour level" ), Qgis::ProcessingNumberParameterType::Double, QVariant(), true
1005 ) );
1006
1007 auto contourLevelList = std::make_unique<QgsProcessingParameterString>(
1008 QStringLiteral( "CONTOUR_LEVEL_LIST" ), QObject::tr( "List of contours level" ), QVariant(), false, true
1009 );
1010 contourLevelList->setHelp( QObject::tr( "Comma separated list of values to export. If filled, the increment, minimum and maximum settings are ignored." ) );
1011 addParameter( contourLevelList.release() );
1012
1013 addParameter( new QgsProcessingParameterCrs( QStringLiteral( "CRS_OUTPUT" ), QObject::tr( "Output coordinate system" ), QVariant(), true ) );
1014
1015
1016 addParameter( new QgsProcessingParameterFeatureSink( QStringLiteral( "OUTPUT_LINES" ), QObject::tr( "Exported contour lines" ), Qgis::ProcessingSourceType::VectorLine ) );
1017 addParameter( new QgsProcessingParameterFeatureSink( QStringLiteral( "OUTPUT_POLYGONS" ), QObject::tr( "Exported contour polygons" ), Qgis::ProcessingSourceType::VectorPolygon ) );
1018}
1019
1020bool QgsMeshContoursAlgorithm::prepareAlgorithm( const QVariantMap &parameters, QgsProcessingContext &context, QgsProcessingFeedback *feedback )
1021{
1022 QgsMeshLayer *meshLayer = parameterAsMeshLayer( parameters, QStringLiteral( "INPUT" ), context );
1023
1024 if ( !meshLayer || !meshLayer->isValid() )
1025 return false;
1026
1027 QgsCoordinateReferenceSystem outputCrs = parameterAsCrs( parameters, QStringLiteral( "CRS_OUTPUT" ), context );
1028 if ( !outputCrs.isValid() )
1029 outputCrs = meshLayer->crs();
1030 mTransform = QgsCoordinateTransform( meshLayer->crs(), outputCrs, context.transformContext() );
1031 if ( !meshLayer->nativeMesh() )
1032 meshLayer->updateTriangularMesh( mTransform ); //necessary to load the native mesh
1033
1034 mTriangularMesh.update( meshLayer->nativeMesh(), mTransform );
1035 mNativeMesh = *meshLayer->nativeMesh();
1036
1037 // Prepare levels
1038 mLevels.clear();
1039 // First, try with the levels list
1040 QString levelsString = parameterAsString( parameters, QStringLiteral( "CONTOUR_LEVEL_LIST" ), context );
1041 if ( !levelsString.isEmpty() )
1042 {
1043 QStringList levelStringList = levelsString.split( ',' );
1044 if ( !levelStringList.isEmpty() )
1045 {
1046 for ( const QString &stringVal : levelStringList )
1047 {
1048 bool ok;
1049 double val = stringVal.toDouble( &ok );
1050 if ( ok )
1051 mLevels.append( val );
1052 else
1053 throw QgsProcessingException( QObject::tr( "Invalid format for level values, must be numbers separated with comma" ) );
1054
1055 if ( mLevels.count() >= 2 )
1056 if ( mLevels.last() <= mLevels.at( mLevels.count() - 2 ) )
1057 throw QgsProcessingException( QObject::tr( "Invalid format for level values, must be different numbers and in increasing order" ) );
1058 }
1059 }
1060 }
1061
1062 if ( mLevels.isEmpty() )
1063 {
1064 double minimum = parameterAsDouble( parameters, QStringLiteral( "MINIMUM" ), context );
1065 double maximum = parameterAsDouble( parameters, QStringLiteral( "MAXIMUM" ), context );
1066 double interval = parameterAsDouble( parameters, QStringLiteral( "INCREMENT" ), context );
1067
1068 if ( interval <= 0 )
1069 throw QgsProcessingException( QObject::tr( "Invalid interval value, must be greater than zero" ) );
1070
1071 if ( minimum >= maximum )
1072 throw QgsProcessingException( QObject::tr( "Invalid minimum and maximum values, minimum must be lesser than maximum" ) );
1073
1074 if ( interval > ( maximum - minimum ) )
1075 throw QgsProcessingException( QObject::tr( "Invalid minimum, maximum and interval values, difference between minimum and maximum must be greater or equal than interval" ) );
1076
1077 int intervalCount = ( maximum - minimum ) / interval;
1078
1079 mLevels.reserve( intervalCount );
1080 for ( int i = 0; i < intervalCount; ++i )
1081 {
1082 mLevels.append( minimum + i * interval );
1083 }
1084 }
1085
1086 // Prepare data
1087 QList<int> datasetGroups = QgsProcessingParameterMeshDatasetGroups::valueAsDatasetGroup( parameters.value( QStringLiteral( "DATASET_GROUPS" ) ) );
1088
1089 if ( feedback )
1090 {
1091 feedback->setProgressText( QObject::tr( "Preparing data" ) );
1092 }
1093
1094 // Extract the date time used to export dataset values under a relative time
1095 QVariant parameterTimeVariant = parameters.value( QStringLiteral( "DATASET_TIME" ) );
1096 QgsInterval relativeTime = datasetRelativetime( parameterTimeVariant, meshLayer, context );
1097
1098 mDateTimeString = meshLayer->formatTime( relativeTime.hours() );
1099
1100 extractDatasetValues( datasetGroups, meshLayer, mNativeMesh, relativeTime, supportedDataType(), mDataPerGroup, feedback );
1101
1102 mLayerRendererSettings = meshLayer->rendererSettings();
1103
1104 return true;
1105}
1106
1107QVariantMap QgsMeshContoursAlgorithm::processAlgorithm( const QVariantMap &parameters, QgsProcessingContext &context, QgsProcessingFeedback *feedback )
1108{
1109 //First, if present, average 3D staked dataset value to 2D face value
1110 const QgsMesh3DAveragingMethod *avgMethod = mLayerRendererSettings.averagingMethod();
1111 for ( DataGroup &dataGroup : mDataPerGroup )
1112 {
1113 if ( dataGroup.dataset3dStakedValue.isValid() )
1114 dataGroup.datasetValues = avgMethod->calculate( dataGroup.dataset3dStakedValue );
1115 }
1116
1117 // Create vector layers
1118 QgsFields polygonFields;
1119 QgsFields lineFields;
1120 polygonFields.append( QgsField( QObject::tr( "group" ), QMetaType::Type::QString ) );
1121 polygonFields.append( QgsField( QObject::tr( "time" ), QMetaType::Type::QString ) );
1122 polygonFields.append( QgsField( QObject::tr( "min_value" ), QMetaType::Type::Double ) );
1123 polygonFields.append( QgsField( QObject::tr( "max_value" ), QMetaType::Type::Double ) );
1124 lineFields.append( QgsField( QObject::tr( "group" ), QMetaType::Type::QString ) );
1125 lineFields.append( QgsField( QObject::tr( "time" ), QMetaType::Type::QString ) );
1126 lineFields.append( QgsField( QObject::tr( "value" ), QMetaType::Type::Double ) );
1127
1128 QgsCoordinateReferenceSystem outputCrs = parameterAsCrs( parameters, QStringLiteral( "CRS_OUTPUT" ), context );
1129
1130 QString lineIdentifier;
1131 QString polygonIdentifier;
1132 std::unique_ptr<QgsFeatureSink> sinkPolygons( parameterAsSink(
1133 parameters,
1134 QStringLiteral( "OUTPUT_POLYGONS" ),
1135 context,
1136 polygonIdentifier,
1137 polygonFields,
1139 outputCrs
1140 ) );
1141 std::unique_ptr<QgsFeatureSink> sinkLines( parameterAsSink(
1142 parameters,
1143 QStringLiteral( "OUTPUT_LINES" ),
1144 context,
1145 lineIdentifier,
1146 lineFields,
1148 outputCrs
1149 ) );
1150
1151 if ( !sinkLines || !sinkPolygons )
1152 return QVariantMap();
1153
1154
1155 for ( int i = 0; i < mDataPerGroup.count(); ++i )
1156 {
1157 DataGroup dataGroup = mDataPerGroup.at( i );
1158 bool scalarDataOnVertices = dataGroup.metadata.dataType() == QgsMeshDatasetGroupMetadata::DataOnVertices;
1159 int count = scalarDataOnVertices ? mNativeMesh.vertices.count() : mNativeMesh.faces.count();
1160
1161 QVector<double> values;
1162 if ( dataGroup.datasetValues.isValid() )
1163 {
1164 // vals could be scalar or vectors, for contour rendering we want always magnitude
1165 values = QgsMeshLayerUtils::calculateMagnitudes( dataGroup.datasetValues );
1166 }
1167 else
1168 {
1169 values = QVector<double>( count, std::numeric_limits<double>::quiet_NaN() );
1170 }
1171
1172 if ( ( !scalarDataOnVertices ) )
1173 {
1174 values = QgsMeshLayerUtils::interpolateFromFacesData(
1175 values,
1176 mNativeMesh,
1177 &dataGroup.activeFaces,
1179 );
1180 }
1181
1182 QgsMeshContours contoursExported( mTriangularMesh, mNativeMesh, values, dataGroup.activeFaces );
1183
1184 QgsAttributes firstAttributes;
1185 firstAttributes.append( dataGroup.metadata.name() );
1186 firstAttributes.append( mDateTimeString );
1187
1188 for ( double level : std::as_const( mLevels ) )
1189 {
1190 QgsGeometry line = contoursExported.exportLines( level, feedback );
1191 if ( feedback->isCanceled() )
1192 return QVariantMap();
1193 if ( line.isEmpty() )
1194 continue;
1195 QgsAttributes lineAttributes = firstAttributes;
1196 lineAttributes.append( level );
1197
1198 QgsFeature lineFeat;
1199 lineFeat.setGeometry( line );
1200 lineFeat.setAttributes( lineAttributes );
1201
1202 if ( !sinkLines->addFeature( lineFeat, QgsFeatureSink::FastInsert ) )
1203 throw QgsProcessingException( writeFeatureError( sinkLines.get(), parameters, QStringLiteral( "OUTPUT_LINES" ) ) );
1204 }
1205
1206 for ( int l = 0; l < mLevels.count() - 1; ++l )
1207 {
1208 QgsGeometry polygon = contoursExported.exportPolygons( mLevels.at( l ), mLevels.at( l + 1 ), feedback );
1209 if ( feedback->isCanceled() )
1210 return QVariantMap();
1211
1212 if ( polygon.isEmpty() )
1213 continue;
1214 QgsAttributes polygonAttributes = firstAttributes;
1215 polygonAttributes.append( mLevels.at( l ) );
1216 polygonAttributes.append( mLevels.at( l + 1 ) );
1217
1218 QgsFeature polygonFeature;
1219 polygonFeature.setGeometry( polygon );
1220 polygonFeature.setAttributes( polygonAttributes );
1221 sinkPolygons->addFeature( polygonFeature );
1222 }
1223
1224 if ( feedback )
1225 {
1226 feedback->setProgress( 100 * i / mDataPerGroup.count() );
1227 }
1228 }
1229
1230 if ( sinkPolygons )
1231 sinkPolygons->finalize();
1232 if ( sinkLines )
1233 sinkLines->finalize();
1234
1235 QVariantMap ret;
1236 ret[QStringLiteral( "OUTPUT_LINES" )] = lineIdentifier;
1237 ret[QStringLiteral( "OUTPUT_POLYGONS" )] = polygonIdentifier;
1238
1239 return ret;
1240}
1241
1242QString QgsMeshExportCrossSection::name() const
1243{
1244 return QStringLiteral( "meshexportcrosssection" );
1245}
1246
1247QString QgsMeshExportCrossSection::displayName() const
1248{
1249 return QObject::tr( "Export cross section dataset values on lines from mesh" );
1250}
1251
1252QString QgsMeshExportCrossSection::group() const
1253{
1254 return QObject::tr( "Mesh" );
1255}
1256
1257QString QgsMeshExportCrossSection::groupId() const
1258{
1259 return QStringLiteral( "mesh" );
1260}
1261
1262QString QgsMeshExportCrossSection::shortHelpString() const
1263{
1264 return QObject::tr( "This algorithm extracts mesh's dataset values from line contained in a vector layer.\n"
1265 "Each line is discretized with a resolution distance parameter for extraction of values on its vertices." );
1266}
1267
1268QString QgsMeshExportCrossSection::shortDescription() const
1269{
1270 return QObject::tr( "Extracts a mesh dataset's values from lines contained in a vector layer." );
1271}
1272
1273QgsProcessingAlgorithm *QgsMeshExportCrossSection::createInstance() const
1274{
1275 return new QgsMeshExportCrossSection();
1276}
1277
1278void QgsMeshExportCrossSection::initAlgorithm( const QVariantMap &configuration )
1279{
1280 Q_UNUSED( configuration );
1281
1282 addParameter( new QgsProcessingParameterMeshLayer( QStringLiteral( "INPUT" ), QObject::tr( "Input mesh layer" ) ) );
1283
1285 QStringLiteral( "DATASET_GROUPS" ),
1286 QObject::tr( "Dataset groups" ),
1287 QStringLiteral( "INPUT" ),
1288 supportedDataType()
1289 ) );
1290
1291 addParameter( new QgsProcessingParameterMeshDatasetTime(
1292 QStringLiteral( "DATASET_TIME" ),
1293 QObject::tr( "Dataset time" ),
1294 QStringLiteral( "INPUT" ),
1295 QStringLiteral( "DATASET_GROUPS" )
1296 ) );
1297
1298 QList<int> datatype;
1299 datatype << static_cast<int>( Qgis::ProcessingSourceType::VectorLine );
1300 addParameter( new QgsProcessingParameterFeatureSource(
1301 QStringLiteral( "INPUT_LINES" ), QObject::tr( "Lines for data export" ), datatype, QVariant(), false
1302 ) );
1303
1304 addParameter( new QgsProcessingParameterDistance(
1305 QStringLiteral( "RESOLUTION" ), QObject::tr( "Line segmentation resolution" ), 10.0, QStringLiteral( "INPUT_LINES" ), false, 0
1306 ) );
1307
1308 addParameter( new QgsProcessingParameterNumber(
1309 QStringLiteral( "COORDINATES_DIGITS" ), QObject::tr( "Digits count for coordinates" ), Qgis::ProcessingNumberParameterType::Integer, 2
1310 ) );
1311
1312 addParameter( new QgsProcessingParameterNumber(
1313 QStringLiteral( "DATASET_DIGITS" ), QObject::tr( "Digits count for dataset value" ), Qgis::ProcessingNumberParameterType::Integer, 2
1314 ) );
1315
1316 addParameter( new QgsProcessingParameterFileDestination(
1317 QStringLiteral( "OUTPUT" ), QObject::tr( "Exported data CSV file" ), QObject::tr( "CSV file (*.csv)" )
1318 ) );
1319}
1320
1321bool QgsMeshExportCrossSection::prepareAlgorithm( const QVariantMap &parameters, QgsProcessingContext &context, QgsProcessingFeedback *feedback )
1322{
1323 QgsMeshLayer *meshLayer = parameterAsMeshLayer( parameters, QStringLiteral( "INPUT" ), context );
1324
1325 if ( !meshLayer || !meshLayer->isValid() )
1326 return false;
1327
1328 mMeshLayerCrs = meshLayer->crs();
1329 mTriangularMesh.update( meshLayer->nativeMesh() );
1330 QList<int> datasetGroups = QgsProcessingParameterMeshDatasetGroups::valueAsDatasetGroup( parameters.value( QStringLiteral( "DATASET_GROUPS" ) ) );
1331
1332 if ( feedback )
1333 {
1334 feedback->setProgressText( QObject::tr( "Preparing data" ) );
1335 }
1336
1337 // Extract the date time used to export dataset values under a relative time
1338 QVariant parameterTimeVariant = parameters.value( QStringLiteral( "DATASET_TIME" ) );
1339 QgsInterval relativeTime = datasetRelativetime( parameterTimeVariant, meshLayer, context );
1340
1341 extractDatasetValues( datasetGroups, meshLayer, *meshLayer->nativeMesh(), relativeTime, supportedDataType(), mDataPerGroup, feedback );
1342
1343 mLayerRendererSettings = meshLayer->rendererSettings();
1344
1345 return true;
1346}
1347
1348QVariantMap QgsMeshExportCrossSection::processAlgorithm( const QVariantMap &parameters, QgsProcessingContext &context, QgsProcessingFeedback *feedback )
1349{
1350 if ( feedback )
1351 feedback->setProgress( 0 );
1352 //First, if present, average 3D staked dataset value to 2D face value
1353 const QgsMesh3DAveragingMethod *avgMethod = mLayerRendererSettings.averagingMethod();
1354 for ( DataGroup &dataGroup : mDataPerGroup )
1355 {
1356 if ( dataGroup.dataset3dStakedValue.isValid() )
1357 dataGroup.datasetValues = avgMethod->calculate( dataGroup.dataset3dStakedValue );
1358 }
1359 double resolution = parameterAsDouble( parameters, QStringLiteral( "RESOLUTION" ), context );
1360 int datasetDigits = parameterAsInt( parameters, QStringLiteral( "DATASET_DIGITS" ), context );
1361 int coordDigits = parameterAsInt( parameters, QStringLiteral( "COORDINATES_DIGITS" ), context );
1362
1363 std::unique_ptr<QgsProcessingFeatureSource> featureSource( parameterAsSource( parameters, QStringLiteral( "INPUT_LINES" ), context ) );
1364 if ( !featureSource )
1365 throw QgsProcessingException( QObject::tr( "Input lines vector layer required" ) );
1366
1367 QgsCoordinateTransform transform( featureSource->sourceCrs(), mMeshLayerCrs, context.transformContext() );
1368
1369 QString outputFileName = parameterAsFileOutput( parameters, QStringLiteral( "OUTPUT" ), context );
1370 QFile file( outputFileName );
1371 if ( !file.open( QIODevice::WriteOnly | QIODevice::Truncate ) )
1372 throw QgsProcessingException( QObject::tr( "Unable to create the output file" ) );
1373
1374 QTextStream textStream( &file );
1375#if QT_VERSION < QT_VERSION_CHECK( 6, 0, 0 )
1376 textStream.setCodec( "UTF-8" );
1377#endif
1378 QStringList header;
1379 header << QStringLiteral( "fid" ) << QStringLiteral( "x" ) << QStringLiteral( "y" ) << QObject::tr( "offset" );
1380 for ( const DataGroup &datagroup : std::as_const( mDataPerGroup ) )
1381 header << datagroup.metadata.name();
1382 textStream << header.join( ',' ) << QStringLiteral( "\n" );
1383
1384 long long featCount = featureSource->featureCount();
1385 long long featCounter = 0;
1386 QgsFeatureIterator featIt = featureSource->getFeatures();
1387 QgsFeature feat;
1388 while ( featIt.nextFeature( feat ) )
1389 {
1390 QgsFeatureId fid = feat.id();
1391 QgsGeometry line = feat.geometry();
1392 try
1393 {
1394 line.transform( transform );
1395 }
1396 catch ( QgsCsException & )
1397 {
1398 line = feat.geometry();
1399 feedback->reportError( QObject::tr( "Could not transform line to mesh CRS" ) );
1400 }
1401
1402 if ( line.isEmpty() )
1403 continue;
1404 double offset = 0;
1405 while ( offset <= line.length() )
1406 {
1407 if ( feedback->isCanceled() )
1408 return QVariantMap();
1409
1410 QStringList textLine;
1411 QgsPointXY point = line.interpolate( offset ).asPoint();
1412 int triangularFaceIndex = mTriangularMesh.faceIndexForPoint_v2( point );
1413 textLine << QString::number( fid ) << QString::number( point.x(), 'f', coordDigits ) << QString::number( point.y(), 'f', coordDigits ) << QString::number( offset, 'f', coordDigits );
1414 if ( triangularFaceIndex >= 0 )
1415 {
1416 //extract dataset values for the point
1417 QgsAttributes attributes;
1418 int nativeFaceIndex = mTriangularMesh.trianglesToNativeFaces().at( triangularFaceIndex );
1419 for ( int i = 0; i < mDataPerGroup.count(); ++i )
1420 {
1421 const DataGroup &dataGroup = mDataPerGroup.at( i );
1422 bool faceActive = dataGroup.activeFaces.active( nativeFaceIndex );
1423 if ( !faceActive )
1424 continue;
1425 QgsMeshDatasetValue value = extractDatasetValue(
1426 point,
1427 nativeFaceIndex,
1428 triangularFaceIndex,
1429 mTriangularMesh,
1430 dataGroup.activeFaces,
1431 dataGroup.datasetValues,
1432 dataGroup.metadata
1433 );
1434
1435 if ( abs( value.x() ) == std::numeric_limits<double>::quiet_NaN() )
1436 textLine << QString( ' ' );
1437 else
1438 textLine << QString::number( value.scalar(), 'f', datasetDigits );
1439 }
1440 }
1441 else
1442 for ( int i = 0; i < mDataPerGroup.count(); ++i )
1443 textLine << QString( ' ' );
1444
1445 textStream << textLine.join( ',' ) << QStringLiteral( "\n" );
1446
1447 offset += resolution;
1448 }
1449
1450 if ( feedback )
1451 {
1452 feedback->setProgress( 100.0 * featCounter / featCount );
1453 if ( feedback->isCanceled() )
1454 return QVariantMap();
1455 }
1456 }
1457
1458 file.close();
1459
1460 QVariantMap ret;
1461 ret[QStringLiteral( "OUTPUT" )] = outputFileName;
1462 return ret;
1463}
1464
1465QString QgsMeshExportTimeSeries::name() const
1466{
1467 return QStringLiteral( "meshexporttimeseries" );
1468}
1469
1470QString QgsMeshExportTimeSeries::displayName() const
1471{
1472 return QObject::tr( "Export time series values from points of a mesh dataset" );
1473}
1474
1475QString QgsMeshExportTimeSeries::group() const
1476{
1477 return QObject::tr( "Mesh" );
1478}
1479
1480QString QgsMeshExportTimeSeries::groupId() const
1481{
1482 return QStringLiteral( "mesh" );
1483}
1484
1485QString QgsMeshExportTimeSeries::shortHelpString() const
1486{
1487 return QObject::tr( "This algorithm extracts mesh's dataset time series values from points contained in a vector layer.\n"
1488 "If the time step is kept to its default value (0 hours), the time step used is the one of the two first datasets of the first selected dataset group." );
1489}
1490
1491QString QgsMeshExportTimeSeries::shortDescription() const
1492{
1493 return QObject::tr( "Extracts a mesh dataset's time series values from points contained in a vector layer." );
1494}
1495
1496QgsProcessingAlgorithm *QgsMeshExportTimeSeries::createInstance() const
1497{
1498 return new QgsMeshExportTimeSeries();
1499}
1500
1501void QgsMeshExportTimeSeries::initAlgorithm( const QVariantMap &configuration )
1502{
1503 Q_UNUSED( configuration );
1504
1505 addParameter( new QgsProcessingParameterMeshLayer( QStringLiteral( "INPUT" ), QObject::tr( "Input mesh layer" ) ) );
1506
1508 QStringLiteral( "DATASET_GROUPS" ),
1509 QObject::tr( "Dataset groups" ),
1510 QStringLiteral( "INPUT" ),
1511 supportedDataType()
1512 ) );
1513
1514 addParameter( new QgsProcessingParameterMeshDatasetTime(
1515 QStringLiteral( "STARTING_TIME" ),
1516 QObject::tr( "Starting time" ),
1517 QStringLiteral( "INPUT" ),
1518 QStringLiteral( "DATASET_GROUPS" )
1519 ) );
1520
1521 addParameter( new QgsProcessingParameterMeshDatasetTime(
1522 QStringLiteral( "FINISHING_TIME" ),
1523 QObject::tr( "Finishing time" ),
1524 QStringLiteral( "INPUT" ),
1525 QStringLiteral( "DATASET_GROUPS" )
1526 ) );
1527
1528 addParameter( new QgsProcessingParameterNumber(
1529 QStringLiteral( "TIME_STEP" ), QObject::tr( "Time step (hours)" ), Qgis::ProcessingNumberParameterType::Double, 0, true, 0
1530 ) );
1531
1532 QList<int> datatype;
1533 datatype << static_cast<int>( Qgis::ProcessingSourceType::VectorPoint );
1534 addParameter( new QgsProcessingParameterFeatureSource(
1535 QStringLiteral( "INPUT_POINTS" ), QObject::tr( "Points for data export" ), datatype, QVariant(), false
1536 ) );
1537
1538 addParameter( new QgsProcessingParameterNumber(
1539 QStringLiteral( "COORDINATES_DIGITS" ), QObject::tr( "Digits count for coordinates" ), Qgis::ProcessingNumberParameterType::Integer, 2
1540 ) );
1541
1542 addParameter( new QgsProcessingParameterNumber(
1543 QStringLiteral( "DATASET_DIGITS" ), QObject::tr( "Digits count for dataset value" ), Qgis::ProcessingNumberParameterType::Integer, 2
1544 ) );
1545
1546 addParameter( new QgsProcessingParameterFileDestination(
1547 QStringLiteral( "OUTPUT" ), QObject::tr( "Exported data CSV file" ), QObject::tr( "CSV file (*.csv)" )
1548 ) );
1549}
1550
1551bool QgsMeshExportTimeSeries::prepareAlgorithm( const QVariantMap &parameters, QgsProcessingContext &context, QgsProcessingFeedback *feedback )
1552{
1553 QgsMeshLayer *meshLayer = parameterAsMeshLayer( parameters, QStringLiteral( "INPUT" ), context );
1554
1555 if ( !meshLayer || !meshLayer->isValid() )
1556 return false;
1557
1558 mMeshLayerCrs = meshLayer->crs();
1559 mTriangularMesh.update( meshLayer->nativeMesh() );
1560
1561 QList<int> datasetGroups = QgsProcessingParameterMeshDatasetGroups::valueAsDatasetGroup( parameters.value( QStringLiteral( "DATASET_GROUPS" ) ) );
1562
1563 if ( feedback )
1564 {
1565 feedback->setProgressText( QObject::tr( "Preparing data" ) );
1566 }
1567
1568 // Extract the date times used to export dataset values
1569 QVariant parameterStartTimeVariant = parameters.value( QStringLiteral( "STARTING_TIME" ) );
1570 QgsInterval relativeStartTime = datasetRelativetime( parameterStartTimeVariant, meshLayer, context );
1571
1572 QVariant parameterEndTimeVariant = parameters.value( QStringLiteral( "FINISHING_TIME" ) );
1573 QgsInterval relativeEndTime = datasetRelativetime( parameterEndTimeVariant, meshLayer, context );
1574
1575 // calculate time steps
1576 qint64 timeStepInterval = parameterAsDouble( parameters, QStringLiteral( "TIME_STEP" ), context ) * 1000 * 3600;
1577 if ( timeStepInterval == 0 )
1578 {
1579 //take the first time step of the first temporal dataset group
1580 for ( int groupIndex : datasetGroups )
1581 {
1582 QgsMeshDatasetGroupMetadata meta = meshLayer->datasetGroupMetadata( QgsMeshDatasetIndex( groupIndex, 0 ) );
1583 if ( !meta.isTemporal() && meshLayer->datasetCount( QgsMeshDatasetIndex( groupIndex, 0 ) ) < 2 )
1584 continue;
1585 else
1586 {
1587 timeStepInterval = meshLayer->datasetRelativeTimeInMilliseconds( QgsMeshDatasetIndex( groupIndex, 1 ) )
1588 - meshLayer->datasetRelativeTimeInMilliseconds( QgsMeshDatasetIndex( groupIndex, 0 ) );
1589 break;
1590 }
1591 }
1592 }
1593
1594 mRelativeTimeSteps.clear();
1595 mTimeStepString.clear();
1596 if ( timeStepInterval != 0 )
1597 {
1598 mRelativeTimeSteps.append( relativeStartTime.seconds() * 1000 );
1599 while ( mRelativeTimeSteps.last() < relativeEndTime.seconds() * 1000 )
1600 mRelativeTimeSteps.append( mRelativeTimeSteps.last() + timeStepInterval );
1601
1602 for ( qint64 relativeTimeStep : std::as_const( mRelativeTimeSteps ) )
1603 {
1604 mTimeStepString.append( meshLayer->formatTime( relativeTimeStep / 3600.0 / 1000.0 ) );
1605 }
1606 }
1607
1608 //Extract needed dataset values
1609 for ( int i = 0; i < datasetGroups.count(); ++i )
1610 {
1611 int groupIndex = datasetGroups.at( i );
1612 QgsMeshDatasetGroupMetadata meta = meshLayer->datasetGroupMetadata( QgsMeshDatasetIndex( groupIndex, 0 ) );
1613 if ( supportedDataType().contains( meta.dataType() ) )
1614 {
1615 mGroupIndexes.append( groupIndex );
1616 mGroupsMetadata[groupIndex] = meta;
1617 int valueCount = meta.dataType() == QgsMeshDatasetGroupMetadata::DataOnVertices ? mTriangularMesh.vertices().count() : meshLayer->nativeMesh()->faceCount();
1618
1619 if ( !mRelativeTimeSteps.isEmpty() )
1620 {
1621 //QMap<qint64, DataGroup> temporalGroup;
1622 QgsMeshDatasetIndex lastDatasetIndex;
1623 for ( qint64 relativeTimeStep : std::as_const( mRelativeTimeSteps ) )
1624 {
1625 QMap<int, int> &groupIndexToData = mRelativeTimeToData[relativeTimeStep];
1626 QgsInterval timeStepInterval( relativeTimeStep / 1000.0 );
1627 QgsMeshDatasetIndex datasetIndex = meshLayer->datasetIndexAtRelativeTime( timeStepInterval, groupIndex );
1628 if ( !datasetIndex.isValid() )
1629 continue;
1630 if ( datasetIndex != lastDatasetIndex )
1631 {
1632 DataGroup dataGroup;
1633 dataGroup.metadata = meta;
1634 dataGroup.datasetValues = meshLayer->datasetValues( datasetIndex, 0, valueCount );
1635 dataGroup.activeFaces = meshLayer->areFacesActive( datasetIndex, 0, meshLayer->nativeMesh()->faceCount() );
1636 if ( dataGroup.metadata.dataType() == QgsMeshDatasetGroupMetadata::DataOnVolumes )
1637 {
1638 dataGroup.dataset3dStakedValue = meshLayer->dataset3dValues( datasetIndex, 0, valueCount );
1639 }
1640 mDatasets.append( dataGroup );
1641 lastDatasetIndex = datasetIndex;
1642 }
1643 groupIndexToData[groupIndex] = mDatasets.count() - 1;
1644 }
1645 }
1646 else
1647 {
1648 // we have only static dataset group
1649 QMap<int, int> &groupIndexToData = mRelativeTimeToData[0];
1650 QgsMeshDatasetIndex datasetIndex( groupIndex, 0 );
1651 DataGroup dataGroup;
1652 dataGroup.metadata = meta;
1653 dataGroup.datasetValues = meshLayer->datasetValues( datasetIndex, 0, valueCount );
1654 dataGroup.activeFaces = meshLayer->areFacesActive( datasetIndex, 0, meshLayer->nativeMesh()->faceCount() );
1655 if ( dataGroup.metadata.dataType() == QgsMeshDatasetGroupMetadata::DataOnVolumes )
1656 {
1657 dataGroup.dataset3dStakedValue = meshLayer->dataset3dValues( datasetIndex, 0, valueCount );
1658 }
1659 mDatasets.append( dataGroup );
1660 groupIndexToData[groupIndex] = mDatasets.count() - 1;
1661 }
1662 }
1663
1664 if ( feedback )
1665 feedback->setProgress( 100 * i / datasetGroups.count() );
1666 }
1667
1668 mLayerRendererSettings = meshLayer->rendererSettings();
1669
1670 return true;
1671}
1672
1673
1674QVariantMap QgsMeshExportTimeSeries::processAlgorithm( const QVariantMap &parameters, QgsProcessingContext &context, QgsProcessingFeedback *feedback )
1675{
1676 if ( feedback )
1677 feedback->setProgress( 0 );
1678 //First, if present, average 3D staked dataset value to 2D face value
1679 const QgsMesh3DAveragingMethod *avgMethod = mLayerRendererSettings.averagingMethod();
1680
1681 for ( DataGroup &dataGroup : mDatasets )
1682 {
1683 if ( dataGroup.dataset3dStakedValue.isValid() )
1684 dataGroup.datasetValues = avgMethod->calculate( dataGroup.dataset3dStakedValue );
1685 }
1686
1687 int datasetDigits = parameterAsInt( parameters, QStringLiteral( "DATASET_DIGITS" ), context );
1688 int coordDigits = parameterAsInt( parameters, QStringLiteral( "COORDINATES_DIGITS" ), context );
1689
1690 std::unique_ptr<QgsProcessingFeatureSource> featureSource( parameterAsSource( parameters, QStringLiteral( "INPUT_POINTS" ), context ) );
1691 if ( !featureSource )
1692 throw QgsProcessingException( QObject::tr( "Input points vector layer required" ) );
1693
1694 QgsCoordinateTransform transform( featureSource->sourceCrs(), mMeshLayerCrs, context.transformContext() );
1695
1696 QString outputFileName = parameterAsFileOutput( parameters, QStringLiteral( "OUTPUT" ), context );
1697 QFile file( outputFileName );
1698 if ( !file.open( QIODevice::WriteOnly | QIODevice::Truncate ) )
1699 throw QgsProcessingException( QObject::tr( "Unable to create the output file" ) );
1700
1701 QTextStream textStream( &file );
1702#if QT_VERSION < QT_VERSION_CHECK( 6, 0, 0 )
1703 textStream.setCodec( "UTF-8" );
1704#endif
1705 QStringList header;
1706 header << QStringLiteral( "fid" ) << QStringLiteral( "x" ) << QStringLiteral( "y" ) << QObject::tr( "time" );
1707
1708 for ( int gi : std::as_const( mGroupIndexes ) )
1709 header << mGroupsMetadata.value( gi ).name();
1710
1711 textStream << header.join( ',' ) << QStringLiteral( "\n" );
1712
1713 long long featCount = featureSource->featureCount();
1714 long long featCounter = 0;
1715 QgsFeatureIterator featIt = featureSource->getFeatures();
1716 QgsFeature feat;
1717 while ( featIt.nextFeature( feat ) )
1718 {
1719 QgsFeatureId fid = feat.id();
1720 QgsGeometry geom = feat.geometry();
1721 try
1722 {
1723 geom.transform( transform );
1724 }
1725 catch ( QgsCsException & )
1726 {
1727 geom = feat.geometry();
1728 feedback->reportError( QObject::tr( "Could not transform line to mesh CRS" ) );
1729 }
1730
1731 if ( geom.isEmpty() )
1732 continue;
1733
1734 QgsPointXY point = geom.asPoint();
1735 int triangularFaceIndex = mTriangularMesh.faceIndexForPoint_v2( point );
1736
1737 if ( triangularFaceIndex >= 0 )
1738 {
1739 int nativeFaceIndex = mTriangularMesh.trianglesToNativeFaces().at( triangularFaceIndex );
1740 if ( !mRelativeTimeSteps.isEmpty() )
1741 {
1742 for ( int timeIndex = 0; timeIndex < mRelativeTimeSteps.count(); ++timeIndex )
1743 {
1744 qint64 timeStep = mRelativeTimeSteps.at( timeIndex );
1745 QStringList textLine;
1746 textLine << QString::number( fid )
1747 << QString::number( point.x(), 'f', coordDigits )
1748 << QString::number( point.y(), 'f', coordDigits )
1749 << mTimeStepString.at( timeIndex );
1750
1751 if ( mRelativeTimeToData.contains( timeStep ) )
1752 {
1753 const QMap<int, int> &groupToData = mRelativeTimeToData.value( timeStep );
1754 for ( int groupIndex : std::as_const( mGroupIndexes ) )
1755 {
1756 if ( !groupToData.contains( groupIndex ) )
1757 continue;
1758 int dataIndex = groupToData.value( groupIndex );
1759 if ( dataIndex < 0 || dataIndex > mDatasets.count() - 1 )
1760 continue;
1761
1762 const DataGroup &dataGroup = mDatasets.at( dataIndex );
1763 QgsMeshDatasetValue value = extractDatasetValue( point, nativeFaceIndex, triangularFaceIndex, mTriangularMesh, dataGroup.activeFaces, dataGroup.datasetValues, dataGroup.metadata );
1764 if ( abs( value.x() ) == std::numeric_limits<double>::quiet_NaN() )
1765 textLine << QString( ' ' );
1766 else
1767 textLine << QString::number( value.scalar(), 'f', datasetDigits );
1768 }
1769 }
1770 textStream << textLine.join( ',' ) << QStringLiteral( "\n" );
1771 }
1772 }
1773 else
1774 {
1775 QStringList textLine;
1776 textLine << QString::number( fid )
1777 << QString::number( point.x(), 'f', coordDigits )
1778 << QString::number( point.y(), 'f', coordDigits )
1779 << QObject::tr( "static dataset" );
1780 const QMap<int, int> &groupToData = mRelativeTimeToData.value( 0 );
1781 for ( int groupIndex : std::as_const( mGroupIndexes ) )
1782 {
1783 if ( !groupToData.contains( groupIndex ) )
1784 continue;
1785 int dataIndex = groupToData.value( groupIndex );
1786 if ( dataIndex < 0 || dataIndex > mDatasets.count() - 1 )
1787 continue;
1788 const DataGroup &dataGroup = mDatasets.at( dataIndex );
1789 QgsMeshDatasetValue value = extractDatasetValue( point, nativeFaceIndex, triangularFaceIndex, mTriangularMesh, dataGroup.activeFaces, dataGroup.datasetValues, dataGroup.metadata );
1790 if ( abs( value.x() ) == std::numeric_limits<double>::quiet_NaN() )
1791 textLine << QString( ' ' );
1792 else
1793 textLine << QString::number( value.scalar(), 'f', datasetDigits );
1794 }
1795 textStream << textLine.join( ',' ) << QStringLiteral( "\n" );
1796 }
1797 }
1798 featCounter++;
1799 if ( feedback )
1800 {
1801 feedback->setProgress( 100.0 * featCounter / featCount );
1802 if ( feedback->isCanceled() )
1803 return QVariantMap();
1804 }
1805 }
1806
1807 file.close();
1808
1809 QVariantMap ret;
1810 ret[QStringLiteral( "OUTPUT" )] = outputFileName;
1811 return ret;
1812}
1813
@ VectorPoint
Vector point layers.
Definition qgis.h:3534
@ VectorPolygon
Vector polygon layers.
Definition qgis.h:3536
@ VectorLine
Vector line layers.
Definition qgis.h:3535
@ Float64
Sixty four bit floating point (double).
Definition qgis.h:381
@ Point
Point.
Definition qgis.h:279
@ LineStringZ
LineStringZ.
Definition qgis.h:296
@ PolygonZ
PolygonZ.
Definition qgis.h:297
@ Hidden
Parameter is hidden and should not be shown to users.
Definition qgis.h:3764
@ Advanced
Parameter is an advanced parameter which should be hidden from users by default.
Definition qgis.h:3763
@ Double
Double/float values.
Definition qgis.h:3804
A vector of attributes.
Represents a coordinate reference system (CRS).
bool isValid() const
Returns whether this CRS is correctly initialized and usable.
Handles coordinate transforms between two coordinate systems.
Custom exception class for Coordinate Reference System related exceptions.
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.
@ 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
QgsFeatureId id
Definition qgsfeature.h:66
void setAttributes(const QgsAttributes &attrs)
Sets the feature's attributes.
QgsGeometry geometry
Definition qgsfeature.h:69
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 canceled()
Internal routines can connect to this signal if they use event loop.
void cancel()
Tells the internal routines that the current operation should be canceled. This should be run by the ...
void setProgress(double progress)
Sets the current progress for the feedback object.
Definition qgsfeedback.h:61
Encapsulate a field in an attribute table or data source.
Definition qgsfield.h:54
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:73
A geometry is the spatial representation of a feature.
double length() const
Returns the planar, 2-dimensional length of geometry.
Qgis::GeometryOperationResult transform(const QgsCoordinateTransform &ct, Qgis::TransformDirection direction=Qgis::TransformDirection::Forward, bool transformZ=false)
Transforms this geometry as described by the coordinate transform ct.
QgsGeometry interpolate(double distance) const
Returns an interpolated point on the geometry at the specified distance.
QgsPointXY asPoint() const
Returns the contents of the geometry as a 2-dimensional point.
bool isEmpty() const
Returns true if the geometry is empty (eg a linestring with no vertices, or a collection with no geom...
A representation of the interval between two datetime values.
Definition qgsinterval.h:47
double seconds() const
Returns the interval duration in seconds.
double hours() const
Returns the interval duration in hours.
Line string geometry type, with support for z-dimension and m-values.
QgsCoordinateReferenceSystem crs
Definition qgsmaplayer.h:87
Abstract class for interpolating 3d stacked mesh data to 2d data.
QgsMeshDataBlock calculate(const QgsMesh3DDataBlock &block3d, QgsFeedback *feedback=nullptr) const
Calculated 2d block values from 3d stacked mesh values.
Exporter of contours lines or polygons from a mesh layer.
A block of integers/doubles from a mesh dataset.
QgsMeshDatasetValue value(int index) const
Returns a value represented by the index For active flag the behavior is undefined.
bool active(int index) const
Returns a value for active flag by the index For scalar and vector 2d the behavior is undefined.
int count() const
Number of items stored in the block.
A collection of dataset group metadata such as whether the data is vector or scalar,...
bool isTemporal() const
Returns whether the dataset group is temporal (contains time-related dataset).
bool isVector() const
Returns whether dataset group has vector data.
DataType dataType() const
Returns whether dataset group data is defined on vertices or faces or volumes.
@ DataOnEdges
Data is defined on edges.
@ DataOnFaces
Data is defined on faces.
@ DataOnVertices
Data is defined on vertices.
@ DataOnVolumes
Data is defined on volumes.
An index that identifies the dataset group (e.g.
bool isValid() const
Returns whether index is valid, ie at least groups is set.
Represents a single mesh dataset value.
double y() const
Returns y value.
double scalar() const
Returns magnitude of vector for vector data or scalar value for scalar data.
double x() const
Returns x value.
Implementation of map layer temporal properties for mesh layers.
QDateTime referenceTime() const
Returns the reference time.
Represents a mesh layer supporting display of data on structured or unstructured meshes.
int datasetCount(const QgsMeshDatasetIndex &index) const
Returns the dataset count in the dataset groups.
QgsMeshRendererSettings rendererSettings() const
Returns renderer settings.
void updateTriangularMesh(const QgsCoordinateTransform &transform=QgsCoordinateTransform())
Gets native mesh and updates (creates if it doesn't exist) the base triangular mesh.
QgsMesh * nativeMesh()
Returns native mesh (nullptr before rendering or calling to updateMesh).
QgsMeshDatasetIndex datasetIndexAtRelativeTime(const QgsInterval &relativeTime, int datasetGroupIndex) const
Returns dataset index from datasets group depending on the relative time from the layer reference tim...
QgsMeshDataBlock datasetValues(const QgsMeshDatasetIndex &index, int valueIndex, int count) const
Returns N vector/scalar values from the index from the dataset.
bool isEditable() const override
Returns true if the layer can be edited.
QgsMeshDataBlock areFacesActive(const QgsMeshDatasetIndex &index, int faceIndex, int count) const
Returns whether the faces are active for particular dataset.
QgsInterval datasetRelativeTime(const QgsMeshDatasetIndex &index)
Returns the relative time of the dataset from the reference time of its group.
QgsMapLayerTemporalProperties * temporalProperties() override
Returns the layer's temporal properties.
qint64 datasetRelativeTimeInMilliseconds(const QgsMeshDatasetIndex &index)
Returns the relative time (in milliseconds) of the dataset from the reference time of its group.
QgsMesh3DDataBlock dataset3dValues(const QgsMeshDatasetIndex &index, int faceIndex, int count) const
Returns N vector/scalar values from the face index from the dataset for 3d stacked meshes.
QString formatTime(double hours)
Returns (date) time in hours formatted to human readable form.
QgsMeshDatasetGroupMetadata datasetGroupMetadata(const QgsMeshDatasetIndex &index) const
Returns the dataset groups metadata.
@ NeighbourAverage
Does a simple average of values defined for all surrounding faces/vertices.
static QgsRasterBlock * exportRasterBlock(const QgsMeshLayer &layer, const QgsMeshDatasetIndex &datasetIndex, const QgsCoordinateReferenceSystem &destinationCrs, const QgsCoordinateTransformContext &transformContext, double mapUnitsPerPixel, const QgsRectangle &extent, QgsRasterBlockFeedback *feedback=nullptr)
Exports mesh layer's dataset values as raster block.
Represents a 2D point.
Definition qgspointxy.h:60
double y
Definition qgspointxy.h:64
double x
Definition qgspointxy.h:63
Point geometry type, with support for z-dimension and m-values.
Definition qgspoint.h:49
Abstract base class for processing algorithms.
Contains information about the context in which a processing algorithm is executed.
QgsDateTimeRange currentTimeRange() const
Returns the current time range to use for temporal operations.
QgsCoordinateTransformContext transformContext() const
Returns the coordinate transform context.
Custom exception class for processing related exceptions.
Base class for providing feedback from a processing algorithm.
virtual void reportError(const QString &error, bool fatalError=false)
Reports that the algorithm encountered an error while executing.
virtual void setProgressText(const QString &text)
Sets a progress report text string.
A coordinate reference system parameter for processing algorithms.
A double numeric parameter for distance values.
An enum based parameter for processing algorithms, allowing for selection from predefined values.
A rectangular map extent parameter for processing algorithms.
A feature sink output for processing algorithms.
An input feature source (such as vector layers) parameter for processing algorithms.
A generic file based destination parameter, for specifying the destination path for a file (non-map l...
A parameter for processing algorithms that need a list of mesh dataset groups.
static QList< int > valueAsDatasetGroup(const QVariant &value)
Returns the value as a list if dataset group indexes.
A parameter for processing algorithms that need a list of mesh dataset index from time parameter.
static QString valueAsTimeType(const QVariant &value)
Returns the dataset value time type as a string : current-context-time : the time is store in the pro...
static QgsMeshDatasetIndex timeValueAsDatasetIndex(const QVariant &value)
Returns the value as a QgsMeshDatasetIndex if the value has "dataset-time-step" type.
static QDateTime timeValueAsDefinedDateTime(const QVariant &value)
Returns the value as a QDateTime if the value has "defined-date-time" type.
A mesh layer parameter for processing algorithms.
A numeric parameter for processing algorithms.
A raster layer destination parameter, for specifying the destination path for a raster layer created ...
Feedback object tailored for raster block reading.
The raster file writer which allows you to save a raster to a new file.
static QString driverForExtension(const QString &extension)
Returns the GDAL driver name for a specified file extension.
A rectangle specified with double values.
double xMinimum
double yMinimum
T begin() const
Returns the beginning of the range.
Definition qgsrange.h:446
A triangular/derived mesh with vertices in map coordinates.
const QVector< QgsMeshFace > & triangles() const
Returns triangles.
const QVector< QgsMeshVertex > & vertices() const
Returns vertices in map coordinate system.
bool qgsDoubleNear(double a, double b, double epsilon=4 *std::numeric_limits< double >::epsilon())
Compare two doubles (but allow some difference).
Definition qgis.h:6607
qint64 QgsFeatureId
64 bit feature ids negative numbers are used for uncommitted/newly added features
QVector< int > QgsMeshFace
List of vertex indexes.
QPair< int, int > QgsMeshEdge
Edge is a straight line seqment between 2 points.
Mesh - vertices, edges and faces.
QVector< QgsMeshVertex > vertices
void clear()
Remove all vertices, edges and faces.
int faceCount() const
Returns number of faces.