QGIS API Documentation  3.20.0-Odense (decaadbb31)
qgsmeshtriangulation.cpp
Go to the documentation of this file.
1 /***************************************************************************
2  qgsmeshtriangulation.cpp
3  -----------------
4  begin : August 9th, 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 "qgsmeshtriangulation.h"
20 #include "qgsvectorlayer.h"
22 #include "qgscurve.h"
23 #include "qgscurvepolygon.h"
24 #include "qgsmultisurface.h"
25 #include "qgsmulticurve.h"
26 #include "qgsfeedback.h"
27 #include "qgslogger.h"
28 
30 {
31  mTriangulation.reset( new QgsDualEdgeTriangulation() );
32 }
33 
34 
36 
37 bool QgsMeshTriangulation::addVertices( QgsFeatureIterator &vertexFeatureIterator, int valueAttribute, const QgsCoordinateTransform &transform, QgsFeedback *feedback, long featureCount )
38 {
39  if ( feedback )
40  feedback->setProgress( 0 );
41 
42  QgsFeature feat;
43  long i = 0;
44  while ( vertexFeatureIterator.nextFeature( feat ) )
45  {
46  if ( feedback )
47  {
48  if ( feedback->isCanceled() )
49  break;
50 
51  feedback->setProgress( 100 * i / featureCount );
52  i++;
53  }
54 
55  addVerticesFromFeature( feat, valueAttribute, transform, feedback );
56  }
57 
58  return true;
59 }
60 
61 bool QgsMeshTriangulation::addBreakLines( QgsFeatureIterator &lineFeatureIterator, int valueAttribute, const QgsCoordinateTransform &transform, QgsFeedback *feedback, long featureCount )
62 {
63  if ( feedback )
64  feedback->setProgress( 0 );
65 
66  QgsFeature feat;
67  long i = 0;
68  while ( lineFeatureIterator.nextFeature( feat ) )
69  {
70  if ( feedback )
71  {
72  if ( feedback->isCanceled() )
73  break;
74 
75  feedback->setProgress( 100 * i / featureCount );
76  i++;
77  }
78 
79  QgsWkbTypes::GeometryType geomType = feat.geometry().type();
80  switch ( geomType )
81  {
83  addVerticesFromFeature( feat, valueAttribute, transform, feedback );
84  break;
87  addBreakLinesFromFeature( feat, valueAttribute, transform, feedback );
88  break;
89  default:
90  break;
91  }
92  }
93 
94  return true;
95 }
96 
98 {
99  return mTriangulation->triangulationToMesh( feedback );
100 }
101 
103 {
104  mCrs = crs;
105 }
106 
107 void QgsMeshTriangulation::addVerticesFromFeature( const QgsFeature &feature, int valueAttribute, const QgsCoordinateTransform &transform, QgsFeedback *feedback )
108 {
109  QgsGeometry geom = feature.geometry();
110  try
111  {
112  geom.transform( transform, QgsCoordinateTransform::ForwardTransform, true );
113  }
114  catch ( QgsCsException &cse )
115  {
116  Q_UNUSED( cse )
117  QgsDebugMsgLevel( QStringLiteral( "Caught CRS exception %1" ).arg( cse.what() ), 4 );
118  return;
119  }
120 
122 
123  double value = 0;
124  if ( valueAttribute >= 0 )
125  value = feature.attribute( valueAttribute ).toDouble();
126 
127  while ( vit != geom.vertices_end() )
128  {
129  if ( feedback && feedback->isCanceled() )
130  break;
131  if ( valueAttribute < 0 )
132  mTriangulation->addPoint( *vit );
133  else
134  {
135  mTriangulation->addPoint( QgsPoint( QgsWkbTypes::PointZ, ( *vit ).x(), ( *vit ).y(), value ) );
136  }
137  ++vit;
138  }
139 }
140 
141 void QgsMeshTriangulation::addBreakLinesFromFeature( const QgsFeature &feature, int valueAttribute, const QgsCoordinateTransform &transform, QgsFeedback *feedback )
142 {
143  double valueOnVertex = 0;
144  if ( valueAttribute >= 0 )
145  valueOnVertex = feature.attribute( valueAttribute ).toDouble();
146 
147  //from QgsTinInterpolator::insertData()
148  std::vector<const QgsCurve *> curves;
149  QgsGeometry geom = feature.geometry();
150  try
151  {
152  geom.transform( transform, QgsCoordinateTransform::ForwardTransform, true );
153  }
154  catch ( QgsCsException &cse )
155  {
156  Q_UNUSED( cse )
157  QgsDebugMsgLevel( QStringLiteral( "Caught CRS exception %1" ).arg( cse.what() ), 4 );
158  return;
159  }
160 
162  {
163  std::vector< const QgsCurvePolygon * > polygons;
164  if ( geom.isMultipart() )
165  {
166  const QgsMultiSurface *ms = qgsgeometry_cast< const QgsMultiSurface * >( geom.constGet() );
167  for ( int i = 0; i < ms->numGeometries(); ++i )
168  {
169  polygons.emplace_back( qgsgeometry_cast< const QgsCurvePolygon * >( ms->geometryN( i ) ) );
170  }
171  }
172  else
173  {
174  polygons.emplace_back( qgsgeometry_cast< const QgsCurvePolygon * >( geom.constGet() ) );
175  }
176 
177  for ( const QgsCurvePolygon *polygon : polygons )
178  {
179  if ( feedback && feedback->isCanceled() )
180  break;
181  if ( !polygon )
182  continue;
183 
184  if ( polygon->exteriorRing() )
185  curves.emplace_back( polygon->exteriorRing() );
186 
187  for ( int i = 0; i < polygon->numInteriorRings(); ++i )
188  {
189  if ( feedback && feedback->isCanceled() )
190  break;
191  curves.emplace_back( polygon->interiorRing( i ) );
192  }
193  }
194  }
195  else
196  {
197  if ( geom.isMultipart() )
198  {
199  const QgsMultiCurve *mc = qgsgeometry_cast< const QgsMultiCurve * >( geom.constGet() );
200  for ( int i = 0; i < mc->numGeometries(); ++i )
201  {
202  if ( feedback && feedback->isCanceled() )
203  break;
204  curves.emplace_back( qgsgeometry_cast< const QgsCurve * >( mc->geometryN( i ) ) );
205  }
206  }
207  else
208  {
209  curves.emplace_back( qgsgeometry_cast< const QgsCurve * >( geom.constGet() ) );
210  }
211  }
212 
213  for ( const QgsCurve *curve : curves )
214  {
215  if ( !curve )
216  continue;
217 
218  if ( feedback && feedback->isCanceled() )
219  break;
220 
221  QgsPointSequence linePoints;
222  curve->points( linePoints );
223  bool hasZ = curve->is3D();
224  if ( valueAttribute >= 0 )
225  for ( int i = 0; i < linePoints.count(); ++i )
226  {
227  if ( feedback && feedback->isCanceled() )
228  break;
229  if ( hasZ )
230  linePoints[i].setZ( valueOnVertex );
231  else
232  {
233  const QgsPoint &point = linePoints.at( i );
234  linePoints[i] = QgsPoint( point.x(), point.y(), valueOnVertex );
235  }
236  }
237 
238  mTriangulation->addLine( linePoints, QgsInterpolator::SourceBreakLines );
239  }
240 }
241 
242 QgsMeshZValueDatasetGroup::QgsMeshZValueDatasetGroup( const QString &datasetGroupName, const QgsMesh &mesh ):
243  QgsMeshDatasetGroup( datasetGroupName, QgsMeshDatasetGroupMetadata::DataOnVertices )
244 {
245  mDataset = std::make_unique< QgsMeshZValueDataset >( mesh );
246 }
247 
249 {
251 }
252 
254 {
255  if ( datasetIndex != 0 )
256  return QgsMeshDatasetMetadata();
257 
258  return mDataset->metadata();
259 }
260 
262 
264 {
265  if ( index != 0 )
266  return nullptr;
267 
268  return mDataset.get();
269 }
270 
271 QDomElement QgsMeshZValueDatasetGroup::writeXml( QDomDocument &doc, const QgsReadWriteContext &context ) const
272 {
273  Q_UNUSED( doc );
274  Q_UNUSED( context );
275  return QDomElement();
276 }
277 
279 {
280  for ( const QgsMeshVertex &vertex : mesh.vertices )
281  {
282  if ( vertex.z() < mZMinimum )
283  mZMinimum = vertex.z();
284  if ( vertex.z() > mZMaximum )
285  mZMaximum = vertex.z();
286  }
287 }
288 
290 {
291  if ( valueIndex < 0 || valueIndex >= mMesh.vertexCount() )
292  return QgsMeshDatasetValue();
293 
294  return QgsMeshDatasetValue( mMesh.vertex( valueIndex ).z() );
295 }
296 
297 QgsMeshDataBlock QgsMeshZValueDataset::datasetValues( bool isScalar, int valueIndex, int count ) const
298 {
299  Q_UNUSED( isScalar )
301  QVector<double> zValues( count );
302  for ( int i = valueIndex; i < valueIndex + count; ++i )
303  zValues[i - valueIndex] = mMesh.vertex( i ).z();
304  ret.setValues( zValues );
305  return ret;
306 }
307 
309 {
310  Q_UNUSED( faceIndex );
311  Q_UNUSED( count );
313  ret.setValid( true );
314  return ret;
315 }
316 
317 bool QgsMeshZValueDataset::isActive( int faceIndex ) const
318 {
319  return ( faceIndex > 0 && faceIndex < mMesh.faceCount() );
320 }
321 
323 {
324  return QgsMeshDatasetMetadata( 0, true, mZMinimum, mZMaximum, 0 );
325 }
326 
328 {
329  return mMesh.vertexCount();
330 }
The vertex_iterator class provides STL-style iterator for vertices.
This class represents a coordinate reference system (CRS).
Class for doing transforms between two map coordinate systems.
@ ForwardTransform
Transform from source to destination CRS.
Custom exception class for Coordinate Reference System related exceptions.
Definition: qgsexception.h:66
Curve polygon geometry type.
Abstract base class for curved geometry type.
Definition: qgscurve.h:36
DualEdgeTriangulation is an implementation of a triangulation class based on the dual edge data struc...
QString what() const
Definition: qgsexception.h:48
Wrapper for iterator of features from vector data provider or vector layer.
bool nextFeature(QgsFeature &f)
The feature class encapsulates a single feature including its unique ID, geometry and a list of field...
Definition: qgsfeature.h:56
QgsGeometry geometry
Definition: qgsfeature.h:67
QVariant attribute(const QString &name) const
Lookup attribute value by attribute name.
Definition: qgsfeature.cpp:302
Base class for feedback objects to be used for cancellation of something running in a worker thread.
Definition: qgsfeedback.h:45
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
int numGeometries() const SIP_HOLDGIL
Returns the number of geometries within the collection.
const QgsAbstractGeometry * geometryN(int n) const
Returns a const reference to a geometry from within the collection.
A geometry is the spatial representation of a feature.
Definition: qgsgeometry.h:124
const QgsAbstractGeometry * constGet() const SIP_HOLDGIL
Returns a non-modifiable (const) reference to the underlying abstract geometry primitive.
QgsWkbTypes::Type wkbType() const SIP_HOLDGIL
Returns type of the geometry as a WKB type (point / linestring / polygon etc.)
bool isMultipart() const SIP_HOLDGIL
Returns true if WKB of the geometry is of WKBMulti* type.
QgsWkbTypes::GeometryType type
Definition: qgsgeometry.h:127
QgsAbstractGeometry::vertex_iterator vertices_begin() const
Returns STL-style iterator pointing to the first vertex of the geometry.
OperationResult transform(const QgsCoordinateTransform &ct, QgsCoordinateTransform::TransformDirection direction=QgsCoordinateTransform::ForwardTransform, bool transformZ=false) SIP_THROW(QgsCsException)
Transforms this geometry as described by the coordinate transform ct.
QgsAbstractGeometry::vertex_iterator vertices_end() const
Returns STL-style iterator pointing to the imaginary vertex after the last vertex of the geometry.
@ SourceBreakLines
Break lines.
QgsMeshDataBlock is a block of integers/doubles that can be used to retrieve: active flags (e....
@ ScalarDouble
Scalar double values.
@ ActiveFlagInteger
Integer boolean flag whether face is active.
void setValues(const QVector< double > &vals)
Sets values.
void setValid(bool valid)
Sets block validity.
QgsMeshDatasetGroupMetadata is a collection of dataset group metadata such as whether the data is vec...
Abstract class that represents a dataset group.
void calculateStatistic()
Calculates the statistics (minimum and maximum)
QgsMeshDatasetMetadata is a collection of mesh dataset metadata such as whether the data is valid or ...
QgsMeshDatasetValue represents single dataset value.
Abstract class that represents a dataset.
bool addBreakLines(QgsFeatureIterator &lineFeatureIterator, int valueAttribute, const QgsCoordinateTransform &transformContext, QgsFeedback *feedback=nullptr, long featureCount=1)
Adds break lines from a vector layer, return true if successful.
void setCrs(const QgsCoordinateReferenceSystem &crs)
Sets the coordinate reference system used for the triangulation.
~QgsMeshTriangulation()
Destructor.
QgsMeshTriangulation()
Constructor.
bool addVertices(QgsFeatureIterator &vertexFeatureIterator, int valueAttribute, const QgsCoordinateTransform &transform, QgsFeedback *feedback=nullptr, long featureCount=1)
Adds vertices to the triangulation from a feature iterator, return true if successful.
QgsMesh triangulatedMesh(QgsFeedback *feedback=nullptr) const
Returns the triangulated mesh.
int datasetCount() const override
Returns the count of datasets in the group.
QDomElement writeXml(QDomDocument &doc, const QgsReadWriteContext &context) const override
Write dataset group information in a DOM element.
void initialize() override
Initialize the dataset group.
QgsMeshDataset * dataset(int index) const override
Returns the dataset with index.
QgsMeshDatasetMetadata datasetMetadata(int datasetIndex) const override
Returns the metadata of the dataset with index datasetIndex.
QgsMeshZValueDatasetGroup(const QString &datasetGroupName, const QgsMesh &mesh)
Constructor.
QgsMeshDataBlock areFacesActive(int faceIndex, int count) const override
Returns whether faces are active.
QgsMeshDatasetMetadata metadata() const override
Returns the metadata of the dataset.
QgsMeshZValueDataset(const QgsMesh &mesh)
Constructor with the mesh.
bool isActive(int faceIndex) const override
Returns whether the face is active.
int valuesCount() const override
Returns the values count.
QgsMeshDatasetValue datasetValue(int valueIndex) const override
Returns the value with index valueIndex.
QgsMeshDataBlock datasetValues(bool isScalar, int valueIndex, int count) const override
Returns count values from valueIndex.
Multi curve geometry collection.
Definition: qgsmulticurve.h:30
Multi surface geometry collection.
Point geometry type, with support for z-dimension and m-values.
Definition: qgspoint.h:49
Q_GADGET double x
Definition: qgspoint.h:52
double z
Definition: qgspoint.h:54
double y
Definition: qgspoint.h:53
The class is used as a container of context for various read/write operations on other objects.
static GeometryType geometryType(Type type) SIP_HOLDGIL
Returns the geometry type for a WKB type, e.g., both MultiPolygon and CurvePolygon would have a Polyg...
Definition: qgswkbtypes.h:938
GeometryType
The geometry types are used to group QgsWkbTypes::Type in a coarse way.
Definition: qgswkbtypes.h:141
QVector< QgsPoint > QgsPointSequence
#define QgsDebugMsgLevel(str, level)
Definition: qgslogger.h:39
const QgsCoordinateReferenceSystem & crs
Mesh - vertices, edges and faces.
int vertexCount() const
Returns number of vertices.
QVector< QgsMeshVertex > vertices
int faceCount() const
Returns number of faces.
QgsMeshVertex vertex(int index) const
Returns a vertex at the index.