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