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