QGIS API Documentation  3.25.0-Master (10b47c2603)
qgstopologicalmesh.cpp
Go to the documentation of this file.
1 /***************************************************************************
2  qgstopologicalmesh.cpp - QgsTopologicalMesh
3 
4  ---------------------
5  begin : 18.6.2021
6  copyright : (C) 2021 by Vincent Cloarec
7  email : vcloarec at gmail dot com
8  ***************************************************************************
9  * *
10  * This program is free software; you can redistribute it and/or modify *
11  * it under the terms of the GNU General Public License as published by *
12  * the Free Software Foundation; either version 2 of the License, or *
13  * (at your option) any later version. *
14  * *
15  ***************************************************************************/
16 
17 #include "qgis.h"
18 #include "qgstopologicalmesh.h"
19 #include "qgsmesheditor.h"
20 #include "qgsmessagelog.h"
21 #include "qgsgeometryutils.h"
22 #include "qgslinestring.h"
23 
24 #include <poly2tri.h>
25 #include <QSet>
26 #include <QQueue>
27 
28 static int vertexPositionInFace( int vertexIndex, const QgsMeshFace &face )
29 {
30  return face.indexOf( vertexIndex );
31 }
32 
33 static int vertexPositionInFace( const QgsMesh &mesh, int vertexIndex, int faceIndex )
34 {
35  if ( faceIndex < 0 || faceIndex >= mesh.faceCount() )
36  return -1;
37 
38  return vertexPositionInFace( vertexIndex, mesh.face( faceIndex ) );
39 }
40 
41 static double crossProduct( int centralVertex, int vertex1, int vertex2, const QgsMesh &mesh )
42 {
43  QgsMeshVertex vc = mesh.vertices.at( centralVertex );
44  QgsMeshVertex v1 = mesh.vertices.at( vertex1 );
45  QgsMeshVertex v2 = mesh.vertices.at( vertex2 );
46 
47  double ux1 = v1.x() - vc.x();
48  double uy1 = v1.y() - vc.y();
49  double vx1 = v2.x() - vc.x();
50  double vy1 = v2.y() - vc.y();
51 
52  return ux1 * vy1 - uy1 * vx1;
53 }
54 
55 
56 QgsMeshVertexCirculator::QgsMeshVertexCirculator( const QgsTopologicalMesh &topologicalMesh, int vertexIndex )
57  : mFaces( topologicalMesh.mMesh->faces )
58  , mFacesNeighborhood( topologicalMesh.mFacesNeighborhood )
59  , mVertexIndex( vertexIndex )
60 {
61  if ( vertexIndex >= 0 && vertexIndex < topologicalMesh.mMesh->vertexCount() )
62  {
63  mCurrentFace = topologicalMesh.mVertexToFace[vertexIndex];
64  mIsValid = vertexPositionInFace( *topologicalMesh.mesh(), vertexIndex, mCurrentFace ) != -1;
65  }
66  else
67  {
68  mIsValid = false;
69  }
70 
71  if ( mIsValid )
72  mLastValidFace = mCurrentFace;
73 }
74 
75 QgsMeshVertexCirculator::QgsMeshVertexCirculator( const QgsTopologicalMesh::TopologicalFaces &topologicalFaces, int faceIndex, int vertexIndex )
76  : mFaces( topologicalFaces.mFaces )
77  , mFacesNeighborhood( topologicalFaces.mFacesNeighborhood )
78  , mVertexIndex( vertexIndex )
79 {
80  const QgsMeshFace &face = topologicalFaces.mFaces.at( faceIndex );
81  mIsValid = vertexPositionInFace( vertexIndex, face ) != -1;
82 
83  mCurrentFace = faceIndex;
84  mLastValidFace = mCurrentFace;
85 }
86 
88  : mFaces( topologicalFaces.mFaces )
89  , mFacesNeighborhood( topologicalFaces.mFacesNeighborhood )
90  , mVertexIndex( vertexIndex )
91 {
92  if ( topologicalFaces.mVerticesToFace.contains( vertexIndex ) )
93  mCurrentFace = topologicalFaces.mVerticesToFace.values( vertexIndex ).first();
94  mLastValidFace = mCurrentFace;
95  mIsValid = mCurrentFace != -1;
96 }
97 
99 {
100  if ( mCurrentFace == -1 )
101  mCurrentFace = mLastValidFace;
102  else
103  {
104  int currentPos = positionInCurrentFace();
105  Q_ASSERT( currentPos != -1 );
106 
107  const QgsMeshFace &currentFace = mFaces.at( mCurrentFace );
108  int faceSize = currentFace.size();
109  mLastValidFace = mCurrentFace;
110  mCurrentFace = mFacesNeighborhood.at( mCurrentFace ).at( ( currentPos + faceSize - 1 ) % currentFace.count() );
111  }
112 
113  return mCurrentFace;
114 }
115 
117 {
118  if ( mCurrentFace == -1 )
119  mCurrentFace = mLastValidFace;
120  else
121  {
122  int currentPos = positionInCurrentFace();
123  Q_ASSERT( currentPos != -1 );
124 
125  const QgsMeshFace &currentFace = mFaces.at( mCurrentFace );
126  int faceSize = currentFace.size();
127  mLastValidFace = mCurrentFace;
128  mCurrentFace = mFacesNeighborhood.at( mCurrentFace ).at( ( currentPos ) % faceSize );
129  }
130 
131  return mCurrentFace;
132 }
133 
135 {
136  return mCurrentFace;
137 }
138 
140 {
141  if ( mCurrentFace != -1 )
142  return mFaces.at( mCurrentFace );
143  else
144  return QgsMeshFace();
145 }
146 
148 {
149  if ( !isValid() )
150  return false;
151 
152  if ( mCurrentFace == -1 )
153  mCurrentFace = mLastValidFace;
154 
155  int firstFace = mCurrentFace;
156 
157  while ( turnClockwise() != -1 && currentFaceIndex() != firstFace ) {}
158  if ( mCurrentFace == firstFace )
159  return false; // a complete turn, so the vertex is not a boundary vertex, something went wrong
160 
161  return ( turnCounterClockwise() != -1 );
162 }
163 
165 {
166  if ( !isValid() )
167  return false;
168 
169  if ( mCurrentFace == -1 )
170  mCurrentFace = mLastValidFace;
171 
172  int firstFace = mCurrentFace;
173 
174  while ( turnCounterClockwise() != -1 && currentFaceIndex() != firstFace ) {}
175  if ( mCurrentFace == firstFace )
176  return false; // a complete turn, so the vertex is not a boundary vertex, something went wrong
177 
178  return ( turnClockwise() != -1 );
179 }
180 
182 {
183  if ( mCurrentFace == -1 )
184  return -1;
185 
186  const QgsMeshFace &face = currentFace();
187 
188  if ( face.isEmpty() )
189  return -1;
190 
191  int vertexPosition = vertexPositionInFace( mVertexIndex, currentFace() );
192 
193  if ( vertexPosition == -1 )
194  return -1;
195 
196  return face.at( ( vertexPosition + 1 ) % face.count() );
197 }
198 
200 {
201  if ( mCurrentFace == -1 )
202  return -1;
203 
204  const QgsMeshFace &face = currentFace();
205 
206  if ( face.isEmpty() )
207  return -1;
208 
209  int vertexPosition = vertexPositionInFace( mVertexIndex, currentFace() );
210 
211  if ( vertexPosition == -1 )
212  return -1;
213 
214  return face.at( ( vertexPosition - 1 + face.count() ) % face.count() );
215 }
216 
218 {
219  return mIsValid;
220 }
221 
223 {
224  QList<int> ret;
225  if ( !isValid() )
226  return ret;
227 
228  if ( mCurrentFace != -1 )
229  ret.append( mCurrentFace );
230  while ( turnCounterClockwise() != ret.first() && currentFaceIndex() != -1 )
231  ret.append( currentFaceIndex() );
232 
233 
234  if ( currentFaceIndex() == -1 ) //we encounter a boundary, restart with other direction
235  {
236  ret.clear();
237  if ( turnClockwise() == -1 )
238  return ret;
239  ret.append( currentFaceIndex() );
240  while ( turnClockwise() != -1 )
241  {
242  ret.append( currentFaceIndex() );
243  }
244  }
245 
246  return ret;
247 }
248 
250 {
251  if ( mDegree != -1 )
252  return mDegree;
253 
254  mDegree = 0;
255  if ( !mIsValid )
256  return mDegree;
257 
258  // if we are on the boundary, we count one more to take account of the circulator will
259  // not cover the last vertex (the other vertex on boundary)
261  mDegree = 2;
262  else
263  mDegree = 1;
264 
265  int firstFace = currentFaceIndex();
266 
267  while ( turnClockwise() != firstFace && currentFaceIndex() != -1 )
268  ++mDegree;
269 
270  return mDegree;
271 }
272 
273 
274 int QgsMeshVertexCirculator::positionInCurrentFace() const
275 {
276  if ( mCurrentFace < 0 || mCurrentFace > mFaces.count() )
277  return -1;
278 
279  return vertexPositionInFace( mVertexIndex, mFaces.at( mCurrentFace ) );
280 }
281 
283 {
284  Changes changes;
285  changes.mFacesToAdd = topologicalFaces.mFaces;
286  changes.mAddedFacesFirstIndex = mMesh->faceCount();
287 
288  changes.mFacesNeighborhoodToAdd.resize( changes.mFacesToAdd.count() );
289  for ( int i = 0; i < changes.mFacesToAdd.count(); ++i )
290  changes.mFacesNeighborhoodToAdd[i] = QVector<int>( changes.mFacesToAdd.at( i ).count(), -1 );
291 
292  for ( int boundary : topologicalFaces.mBoundaries )
293  {
294  //if the boundary is a free vertex in the destination mesh, no need to check
295  if ( mVertexToFace.at( boundary ) == -1 )
296  continue;
297 
298  const QList<int> &linkedFaces = topologicalFaces.mVerticesToFace.values( boundary );
299  for ( int linkedFace : linkedFaces )
300  {
301  QgsMeshVertexCirculator newFacesCirculator( topologicalFaces, linkedFace, boundary );
302  //search for face boundary on clockwise side of new faces
303  newFacesCirculator.goBoundaryClockwise();
304  int oppositeVertexForNewFace = newFacesCirculator.oppositeVertexClockwise();
305  if ( mVertexToFace.at( oppositeVertexForNewFace ) == -1 )
306  continue;
307 
308  QgsMeshVertexCirculator meshCirculator = vertexCirculator( boundary );
309  meshCirculator.goBoundaryCounterClockwise();
310  int oppositeVertexForMeshFace = meshCirculator.oppositeVertexCounterClockwise();
311 
312  const QgsMeshFace &newFaceBoundary = newFacesCirculator.currentFace();
313  int boundaryPositionInNewFace = vertexPositionInFace( boundary, newFaceBoundary );
314 
315  if ( oppositeVertexForMeshFace != oppositeVertexForNewFace )
316  {
317  changes.mFacesNeighborhoodToAdd[newFacesCirculator.currentFaceIndex()][boundaryPositionInNewFace] = -1 ;
318  }
319  else
320  {
321  const QgsMeshFace &meshFaceBoundary = meshCirculator.currentFace();
322  int boundaryPositionInMeshFace = vertexPositionInFace( meshCirculator.oppositeVertexCounterClockwise(), meshFaceBoundary );
323 
324  changes.mNeighborhoodChanges.append( std::array<int, 4>(
325  {
326  meshCirculator.currentFaceIndex(),
327  boundaryPositionInMeshFace,
328  -1,
329  changes.addedFaceIndexInMesh( newFacesCirculator.currentFaceIndex() )
330  } ) );
331 
332  changes.mFacesNeighborhoodToAdd[newFacesCirculator.currentFaceIndex()][boundaryPositionInNewFace] = meshCirculator.currentFaceIndex();
333  }
334  }
335  }
336 
337  for ( int f = 0; f < changes.mFacesToAdd.count(); ++f )
338  for ( int n = 0; n < changes.mFacesToAdd.at( f ).count(); ++n )
339  if ( changes.mFacesNeighborhoodToAdd.at( f ).at( n ) == -1 )
340  changes.mFacesNeighborhoodToAdd[f][n] = changes.addedFaceIndexInMesh( topologicalFaces.mFacesNeighborhood.at( f ).at( n ) );
341 
342  const QList<int> &verticesToFaceToChange = topologicalFaces.mVerticesToFace.keys();
343  for ( const int vtc : verticesToFaceToChange )
344  if ( mVertexToFace.at( vtc ) == -1 )
345  changes.mVerticesToFaceChanges.append( {vtc,
346  mVertexToFace.at( vtc ),
347  changes.addedFaceIndexInMesh( topologicalFaces.mVerticesToFace.values( vtc ).first() ) } );
348 
349  applyChanges( changes );
350 
351  return changes;
352 }
353 
355 {
356  int initialVerticesCount = mMesh->vertices.count();
357  if ( !changes.mVerticesToAdd.empty() )
358  {
359  int newSize = mMesh->vertices.count() + changes.mVerticesToAdd.count();
360  mMesh->vertices.resize( newSize );
361  mVertexToFace.resize( newSize );
362  }
363 
364  if ( !changes.mFacesToAdd.empty() )
365  {
366  int newSize = mMesh->faceCount() + changes.mFacesToAdd.count();
367  mMesh->faces.resize( newSize );
368  mFacesNeighborhood.resize( newSize );
369  }
370 
371  for ( int i = 0; i < changes.mFacesToRemove.count(); ++i )
372  {
373  mMesh->faces[changes.removedFaceIndexInMesh( i )] = QgsMeshFace();
374  mFacesNeighborhood[changes.removedFaceIndexInMesh( i )] = FaceNeighbors();//changes.facesNeighborhoodToRemove[i];
375  }
376 
377  for ( int i = 0; i < changes.mVerticesToRemoveIndexes.count(); ++i )
378  {
379  int vertexIndex = changes.mVerticesToRemoveIndexes.at( i );
380  if ( mVertexToFace.at( vertexIndex ) == -1 )
381  dereferenceAsFreeVertex( vertexIndex );
382  mMesh->vertices[vertexIndex] = QgsMeshVertex();
383  mVertexToFace[vertexIndex] = -1;
384  }
385 
386  for ( int i = 0; i < changes.mVerticesToAdd.count(); ++i )
387  {
388  mMesh->vertices[initialVerticesCount + i] = changes.mVerticesToAdd.at( i );
389  mVertexToFace[initialVerticesCount + i] = changes.mVertexToFaceToAdd.at( i );
390  if ( changes.mVertexToFaceToAdd.at( i ) == -1 )
391  referenceAsFreeVertex( initialVerticesCount + i );
392  }
393 
394  for ( int i = 0; i < changes.mFacesToAdd.count(); ++i )
395  {
396  mMesh->faces[changes.addedFaceIndexInMesh( i )] = changes.mFacesToAdd.at( i );
397  mFacesNeighborhood[changes.addedFaceIndexInMesh( i )] = changes.mFacesNeighborhoodToAdd.at( i );
398  }
399 
400  for ( const std::array<int, 4> neigborChange : std::as_const( changes.mNeighborhoodChanges ) )
401  {
402  const int faceIndex = neigborChange.at( 0 );
403  const int positionInFace = neigborChange.at( 1 );
404  const int valueToApply = neigborChange.at( 3 );
405  mFacesNeighborhood[faceIndex][positionInFace] = valueToApply;
406  }
407 
408  for ( const std::array<int, 3> vertexToFaceChange : std::as_const( changes.mVerticesToFaceChanges ) )
409  {
410  int vertexIndex = vertexToFaceChange.at( 0 );
411  mVertexToFace[vertexToFaceChange.at( 0 )] = vertexToFaceChange.at( 2 );
412 
413  if ( vertexToFaceChange.at( 2 ) == -1 &&
414  vertexToFaceChange.at( 1 ) != -1 &&
415  !mMesh->vertices.at( vertexIndex ).isEmpty() )
416  referenceAsFreeVertex( vertexIndex );
417 
418  if ( vertexToFaceChange.at( 1 ) == -1 && vertexToFaceChange.at( 2 ) != -1 )
419  dereferenceAsFreeVertex( vertexIndex );
420  }
421 
422  for ( int i = 0; i < changes.mChangeCoordinateVerticesIndexes.count(); ++i )
423  {
424  int vertexIndex = changes.mChangeCoordinateVerticesIndexes.at( i );
425  if ( !changes.mNewZValues.isEmpty() )
426  mMesh->vertices[vertexIndex].setZ( changes.mNewZValues.at( i ) );
427  if ( !changes.mNewXYValues.isEmpty() )
428  {
429  const QgsPointXY &pt = changes.mNewXYValues.at( i );
430  mMesh->vertices[vertexIndex].setX( pt.x() );
431  mMesh->vertices[vertexIndex].setY( pt.y() );
432  }
433  }
434 }
435 
437 {
438  for ( const std::array<int, 4> neigborChange : std::as_const( changes.mNeighborhoodChanges ) )
439  {
440  const int faceIndex = neigborChange.at( 0 );
441  const int positionInFace = neigborChange.at( 1 );
442  const int valueToApply = neigborChange.at( 2 );
443  mFacesNeighborhood[faceIndex][positionInFace] = valueToApply;
444  }
445 
446  for ( int i = 0; i < changes.mFacesToRemove.count(); ++i )
447  {
448  mMesh->faces[changes.removedFaceIndexInMesh( i )] = changes.mFacesToRemove.at( i );
449  mFacesNeighborhood[changes.removedFaceIndexInMesh( i )] = changes.mFacesNeighborhoodToRemove.at( i );
450  }
451 
452  for ( int i = 0; i < changes.mVerticesToRemoveIndexes.count(); ++i )
453  {
454  int vertexIndex = changes.mVerticesToRemoveIndexes.at( i );
455  mMesh->vertices[vertexIndex] = changes.mRemovedVertices.at( i );
456  mVertexToFace[vertexIndex] = changes.mVerticesToFaceRemoved.at( i );
457  if ( mVertexToFace.at( vertexIndex ) == -1 )
458  referenceAsFreeVertex( vertexIndex );
459  }
460 
461  int verticesToFaceChangesCount = changes.mVerticesToFaceChanges.count();
462  for ( int i = 0; i < verticesToFaceChangesCount; ++i )
463  {
464  const std::array<int, 3> vertexToFaceChange = changes.mVerticesToFaceChanges.at( verticesToFaceChangesCount - i - 1 );
465  int vertexIndex = vertexToFaceChange.at( 0 );
466  mVertexToFace[vertexIndex] = vertexToFaceChange.at( 1 );
467 
468  if ( vertexToFaceChange.at( 2 ) == -1 && vertexToFaceChange.at( 1 ) != -1 )
469  dereferenceAsFreeVertex( vertexIndex );
470 
471  if ( vertexToFaceChange.at( 1 ) == -1 &&
472  vertexToFaceChange.at( 2 ) != -1 &&
473  !mMesh->vertex( vertexIndex ).isEmpty() )
474  referenceAsFreeVertex( vertexIndex );
475  }
476 
477  if ( !changes.mFacesToAdd.empty() )
478  {
479  int newSize = mMesh->faceCount() - changes.mFacesToAdd.count();
480  mMesh->faces.resize( newSize );
481  mFacesNeighborhood.resize( newSize );
482  }
483 
484  if ( !changes.mVerticesToAdd.isEmpty() )
485  {
486  int newSize = mMesh->vertexCount() - changes.mVerticesToAdd.count();
487 
488  for ( int i = newSize; i < mMesh->vertexCount(); ++i )
489  if ( mVertexToFace.at( i ) == -1 )
490  dereferenceAsFreeVertex( i );
491 
492  mMesh->vertices.resize( newSize );
493  mVertexToFace.resize( newSize );
494  }
495 
496  for ( int i = 0; i < changes.mChangeCoordinateVerticesIndexes.count(); ++i )
497  {
498  int vertexIndex = changes.mChangeCoordinateVerticesIndexes.at( i );
499  if ( !changes.mOldZValues.isEmpty() )
500  mMesh->vertices[vertexIndex].setZ( changes.mOldZValues.at( i ) );
501  if ( !changes.mOldXYValues.isEmpty() )
502  {
503  const QgsPointXY &pt = changes.mOldXYValues.at( i );
504  mMesh->vertices[vertexIndex].setX( pt.x() );
505  mMesh->vertices[vertexIndex].setY( pt.y() );
506  }
507  }
508 }
509 
511 {
512  return QgsMeshVertexCirculator( *this, vertexIndex );
513 }
514 
515 QSet<int> QgsTopologicalMesh::concernedFacesBy( const QList<int> faceIndexes ) const
516 {
517  QSet<int> faces;
518  for ( const int faceIndex : faceIndexes )
519  {
520  const QgsMeshFace &face = mMesh->face( faceIndex );
521  for ( int i = 0; i < face.count(); ++i )
522  faces.unite( qgis::listToSet( facesAroundVertex( face.at( i ) ) ) );
523  }
524  return faces;
525 }
526 
527 void QgsTopologicalMesh::dereferenceAsFreeVertex( int vertexIndex )
528 {
529  mFreeVertices.remove( vertexIndex );
530 }
531 
532 void QgsTopologicalMesh::referenceAsFreeVertex( int vertexIndex )
533 {
534  // QSet used to retrieve free vertex without going through all the vertices container.
535  // But maybe later we could use more sophisticated reference (spatial index?), to retrieve free vertex in an extent
536  mFreeVertices.insert( vertexIndex );
537 }
538 
540 {
541  for ( int faceIndex = 0 ; faceIndex < mMesh->faces.count( ); ++faceIndex )
542  {
543  const QgsMeshFace &face = mMesh->faces.at( faceIndex );
544  const FaceNeighbors &neighborhood = mFacesNeighborhood.at( faceIndex );
545  if ( face.count() != neighborhood.count() )
547  for ( int i = 0; i < face.count(); ++i )
548  {
549  int vertexIndex = face.at( i );
550  // check if each vertices is linked to a face (not free vertex)
551  if ( mVertexToFace.at( vertexIndex ) == -1 )
553 
554  int neighborIndex = neighborhood.at( i );
555  if ( neighborIndex != -1 )
556  {
557  const QgsMeshFace &neighborFace = mMesh->faces.at( neighborIndex );
558  if ( neighborFace.isEmpty() )
560  int neighborSize = neighborFace.size();
561  const FaceNeighbors &neighborhoodOfNeighbor = mFacesNeighborhood.at( neighborIndex );
562  int posInNeighbor = vertexPositionInFace( *mMesh, vertexIndex, neighborIndex );
563  if ( neighborhoodOfNeighbor.isEmpty() || neighborhoodOfNeighbor.at( ( posInNeighbor + neighborSize - 1 ) % neighborSize ) != faceIndex )
565  }
566  }
567  }
568 
569  for ( int vertexIndex = 0; vertexIndex < mMesh->vertexCount(); ++vertexIndex )
570  {
571  if ( mVertexToFace.at( vertexIndex ) != -1 )
572  {
573  if ( !mMesh->face( mVertexToFace.at( vertexIndex ) ).contains( vertexIndex ) )
575 
576  if ( facesAroundVertex( vertexIndex ).count() == 0 )
578  }
579  }
580 
581  return QgsMeshEditingError();
582 }
583 
584 QgsMeshEditingError QgsTopologicalMesh::checkTopology( const QgsMesh &mesh, int maxVerticesPerFace )
585 {
586  QgsMesh temp = mesh;
587  QgsMeshEditingError error;
588  createTopologicalMesh( &temp, maxVerticesPerFace, error );
589  return error;
590 }
591 
593 {
594  return mMesh;
595 }
596 
597 int QgsTopologicalMesh::firstFaceLinked( int vertexIndex ) const
598 {
599  if ( vertexIndex < 0 || vertexIndex >= mMesh->vertexCount() )
600  return -1;
601  return mVertexToFace.at( vertexIndex );
602 }
603 
604 bool QgsTopologicalMesh::isVertexOnBoundary( int vertexIndex ) const
605 {
606  QgsMeshVertexCirculator circulator = vertexCirculator( vertexIndex );
607 
608  if ( circulator.isValid() )
609  return circulator.goBoundaryClockwise();
610 
611  return false;
612 }
613 
614 bool QgsTopologicalMesh::isVertexFree( int vertexIndex ) const
615 {
616  if ( vertexIndex < 0 || vertexIndex >= mMesh->vertexCount() )
617  return false;
618 
619  if ( mMesh->vertices.at( vertexIndex ).isEmpty() )
620  return false;
621 
622  return mVertexToFace.at( vertexIndex ) == -1;
623 }
624 
626 {
627 #if QT_VERSION < QT_VERSION_CHECK(5, 14, 0)
628  return mFreeVertices.values();
629 #else
630  return QList<int>( mFreeVertices.begin(), mFreeVertices.end() );
631 #endif
632 }
633 
635 {
636  // First check if the face is convex and put it counter clockwise
637  // If the index are not well ordered (edges intersect), invalid face --> return false
638  int faceSize = face.count();
639  if ( faceSize < 3 )
641 
642  int direction = 0;
643  for ( int i = 0; i < faceSize; ++i )
644  {
645  int iv0 = face[i];
646  int iv1 = face[( i + 1 ) % faceSize];
647  int iv2 = face[( i + 2 ) % faceSize];
648 
649  if ( iv0 < 0 || iv0 >= mesh->vertexCount() )
651 
652  if ( iv1 < 0 || iv1 >= mesh->vertexCount() )
654 
655  if ( iv2 < 0 || iv2 >= mesh->vertexCount() )
657 
658  const QgsMeshVertex &v0 = mesh->vertices.at( iv0 ) ;
659  const QgsMeshVertex &v1 = mesh->vertices.at( iv1 ) ;
660  const QgsMeshVertex &v2 = mesh->vertices.at( iv2 ) ;
661 
662  if ( v0.isEmpty() )
664 
665  if ( v1.isEmpty() )
667 
668  if ( v2.isEmpty() )
670 
671  double crossProd = crossProduct( iv1, iv0, iv2, *mesh ); //if cross product>0, we have two edges clockwise
672  if ( direction != 0 && crossProd * direction < 0 ) // We have a convex face or a (partially) flat face
674  else if ( crossProd == 0 )
676  else if ( direction == 0 && crossProd != 0 )
677  direction = crossProd / std::fabs( crossProd );
678  }
679 
680  if ( direction > 0 )// clockwise --> reverse the order of the index;
681  {
682  for ( int i = 0; i < faceSize / 2; ++i )
683  {
684  int temp = face[i];
685  face[i] = face.at( faceSize - i - 1 );
686  face[faceSize - i - 1] = temp;
687  }
688  }
689 
691 }
692 
694 {
695  QVector<int> oldToNewIndex( mMesh->vertices.count(), -1 );
696  int verticesTotalCount = mMesh->vertexCount();
697  int oldIndex = 0;
698  int newIndex = 0;
699  while ( oldIndex < verticesTotalCount )
700  {
701  if ( mMesh->vertex( oldIndex ).isEmpty() )
702  {
703  oldIndex++;
704  }
705  else
706  {
707  oldToNewIndex[oldIndex] = newIndex;
708  if ( oldIndex != newIndex )
709  mMesh->vertices[newIndex] = mMesh->vertices[oldIndex];
710  oldToNewIndex[oldIndex] = newIndex;
711  oldIndex++;
712  newIndex++;
713  }
714  }
715 
716  mMesh->vertices.resize( newIndex );
717 
718  oldIndex = 0;
719  newIndex = 0;
720  int facesTotalCount = mMesh->faceCount();
721  while ( oldIndex < facesTotalCount )
722  {
723  if ( mMesh->face( oldIndex ).isEmpty() )
724  oldIndex++;
725  else
726  {
727  if ( oldIndex != newIndex )
728  mMesh->faces[newIndex] = mMesh->faces[oldIndex];
729  QgsMeshFace &face = mMesh->faces[newIndex];
730  for ( int i = 0; i < face.count(); ++i )
731  face[i] = oldToNewIndex[face.at( i )];
732  newIndex++;
733  oldIndex++;
734  }
735  }
736 
737  mMesh->faces.resize( newIndex );
738 
739  mVertexToFace.clear();
740  mFacesNeighborhood.clear();
741 }
742 
744 {
745  QVector<int> oldToNewVerticesIndexes;
746  if ( !renumberVertices( oldToNewVerticesIndexes ) )
747  return false;
748 
749 
750  QVector<int> oldToNewFacesIndexes;
751 
752  if ( !renumberFaces( oldToNewFacesIndexes ) )
753  return false;
754 
755  // If we are here, we can apply the renumbering
756 
757  QVector<QgsMeshVertex> tempVertices( mMesh->vertices.count() );
758  for ( int i = 0; i < oldToNewVerticesIndexes.count(); ++i )
759  {
760  tempVertices[oldToNewVerticesIndexes.at( i )] = mMesh->vertex( i );
761  }
762  mMesh->vertices = tempVertices;
763 
764  QVector<QgsMeshFace> tempFaces( mMesh->faces.count() );
765  for ( int i = 0; i < oldToNewFacesIndexes.count(); ++i )
766  {
767  tempFaces[oldToNewFacesIndexes.at( i )] = mMesh->face( i );
768 
769  QgsMeshFace &face = tempFaces[oldToNewFacesIndexes.at( i )];
770 
771  for ( int fi = 0; fi < face.count(); ++fi )
772  {
773  face[fi] = oldToNewVerticesIndexes.at( face.at( fi ) );
774  }
775  }
776 
777  mMesh->faces = tempFaces;
778 
779  return true;
780 
781 }
782 
783 bool QgsTopologicalMesh::renumberVertices( QVector<int> &oldToNewIndex ) const
784 {
785  std::vector<QgsMeshVertexCirculator> circulators;
786  circulators.reserve( mMesh->vertices.count() );
787  int minDegree = std::numeric_limits<int>::max();
788  int minDegreeVertex = -1;
789 
790  QQueue<int> queue;
791  QSet<int> nonThreadedVertex;
792  oldToNewIndex = QVector<int> ( mMesh->vertexCount(), -1 );
793  for ( int i = 0; i < mMesh->vertexCount(); ++i )
794  {
795  circulators.emplace_back( *this, i );
796  const QgsMeshVertexCirculator &circulator = circulators.back();
797  if ( circulators.back().degree() < minDegree )
798  {
799  minDegreeVertex = circulators.size() - 1;
800  minDegree = circulator.degree();
801  }
802  nonThreadedVertex.insert( i );
803  }
804 
805  auto sortedNeighbor = [ = ]( QList<int> &neighbors, int index )
806  {
807  const QgsMeshVertexCirculator &circ = circulators.at( index );
808 
809  if ( !circ.isValid() )
810  return;
811 
813  neighbors.append( circ.oppositeVertexCounterClockwise() );
814 
815  int firsrFace = circ.currentFaceIndex();
816  do
817  {
818  int neighborIndex = circ.oppositeVertexClockwise();
819  int degree = circulators.at( neighborIndex ).degree();
820  QList<int>::Iterator it = neighbors.begin();
821  while ( it != neighbors.end() )
822  {
823  if ( degree <= circulators.at( *it ).degree() )
824  {
825  neighbors.insert( it, neighborIndex );
826  break;
827  }
828  it++;
829  }
830  if ( it == neighbors.end() )
831  neighbors.append( neighborIndex );
832  }
833  while ( circ.turnClockwise() != firsrFace && circ.currentFaceIndex() != -1 );
834  };
835 
836  int newIndex = 0;
837  int currentVertex = minDegreeVertex;
838  nonThreadedVertex.remove( minDegreeVertex );
839 
840  while ( newIndex < mMesh->vertexCount() )
841  {
842  if ( oldToNewIndex[currentVertex] == -1 )
843  oldToNewIndex[currentVertex] = newIndex++;
844 
845  if ( circulators.at( currentVertex ).isValid() )
846  {
847  QList<int> neighbors;
848  sortedNeighbor( neighbors, currentVertex );
849 
850  for ( const int i : std::as_const( neighbors ) )
851  if ( oldToNewIndex.at( i ) == -1 && nonThreadedVertex.contains( i ) )
852  {
853  queue.enqueue( i );
854  nonThreadedVertex.remove( i );
855  }
856  }
857 
858  if ( queue.isEmpty() )
859  {
860  if ( nonThreadedVertex.isEmpty() && newIndex < mMesh->vertexCount() )
861  return false;
862 
863  const QList<int> remainingVertex = qgis::setToList( nonThreadedVertex );
864  int minRemainingDegree = std::numeric_limits<int>::max();
865  int minRemainingVertex = -1;
866  for ( const int i : remainingVertex )
867  {
868  int degree = circulators.at( i ).degree();
869  if ( degree < minRemainingDegree )
870  {
871  minRemainingDegree = degree;
872  minRemainingVertex = i;
873  }
874  }
875  currentVertex = minRemainingVertex;
876  nonThreadedVertex.remove( currentVertex );
877  }
878  else
879  {
880  currentVertex = queue.dequeue();
881  }
882  }
883 
884  return true;
885 }
886 
887 bool QgsTopologicalMesh::renumberFaces( QVector<int> &oldToNewIndex ) const
888 {
889  QQueue<int> queue;
890  QSet<int> nonThreadedFaces;
891 
892  oldToNewIndex = QVector<int>( mMesh->faceCount(), -1 );
893 
894  QVector<int> faceDegrees( mMesh->faceCount(), 0 );
895 
896  int minDegree = std::numeric_limits<int>::max();
897  int minDegreeFace = -1;
898  for ( int faceIndex = 0; faceIndex < mMesh->faceCount(); ++faceIndex )
899  {
900  const FaceNeighbors &neighbors = mFacesNeighborhood.at( faceIndex );
901 
902  int degree = 0;
903  for ( int n = 0; n < neighbors.size(); ++n )
904  {
905  if ( neighbors.at( n ) != -1 )
906  degree++;
907  }
908 
909  if ( degree < minDegree )
910  {
911  minDegree = degree;
912  minDegreeFace = faceIndex;
913  }
914 
915  faceDegrees[faceIndex] = degree;
916  nonThreadedFaces.insert( faceIndex );
917  }
918 
919  int newIndex = 0;
920  int currentFace = minDegreeFace;
921  nonThreadedFaces.remove( minDegreeFace );
922 
923  auto sortedNeighbor = [ = ]( QList<int> &neighbors, int index )
924  {
925  const FaceNeighbors &neighborhood = mFacesNeighborhood.at( index );
926 
927  for ( int i = 0; i < neighborhood.count(); ++i )
928  {
929  int neighborIndex = neighborhood.at( i );
930  if ( neighborIndex == -1 )
931  continue;
932 
933  int degree = faceDegrees.at( neighborIndex );
934  if ( neighbors.isEmpty() )
935  neighbors.append( neighborIndex );
936  else
937  {
938  QList<int>::Iterator it = neighbors.begin();
939  while ( it != neighbors.end() )
940  {
941  if ( degree <= faceDegrees.at( *it ) )
942  {
943  neighbors.insert( it, neighborIndex );
944  break;
945  }
946  it++;
947  }
948  if ( it == neighbors.end() )
949  neighbors.append( neighborIndex );
950  }
951  }
952  };
953 
954  while ( newIndex < mMesh->faceCount() )
955  {
956  if ( oldToNewIndex[currentFace] == -1 )
957  oldToNewIndex[currentFace] = newIndex++;
958 
959  QList<int> neighbors;
960  sortedNeighbor( neighbors, currentFace );
961 
962  for ( const int i : std::as_const( neighbors ) )
963  if ( oldToNewIndex.at( i ) == -1 && nonThreadedFaces.contains( i ) )
964  {
965  queue.enqueue( i );
966  nonThreadedFaces.remove( i );
967  }
968 
969  if ( queue.isEmpty() )
970  {
971  if ( nonThreadedFaces.isEmpty() && newIndex < mMesh->faceCount() )
972  return false;
973 
974  const QList<int> remainingFace = qgis::setToList( nonThreadedFaces );
975  int minRemainingDegree = std::numeric_limits<int>::max();
976  int minRemainingFace = -1;
977  for ( const int i : remainingFace )
978  {
979  int degree = faceDegrees.at( i );
980  if ( degree < minRemainingDegree )
981  {
982  minRemainingDegree = degree;
983  minRemainingFace = i;
984  }
985  }
986  currentFace = minRemainingFace;
987  nonThreadedFaces.remove( currentFace );
988  }
989  else
990  {
991  currentFace = queue.dequeue();
992  }
993 
994  }
995 
996  return true;
997 }
998 
999 QVector<QgsMeshFace> QgsTopologicalMesh::Changes::addedFaces() const
1000 {
1001  return mFacesToAdd;
1002 }
1003 
1004 QVector<QgsMeshFace> QgsTopologicalMesh::Changes::removedFaces() const
1005 {
1006  return mFacesToRemove;
1007 }
1008 
1010 {
1011  return mFaceIndexesToRemove;
1012 }
1013 
1014 QVector<QgsMeshVertex> QgsTopologicalMesh::Changes::addedVertices() const
1015 {
1016  return mVerticesToAdd;
1017 }
1018 
1020 {
1021  return mChangeCoordinateVerticesIndexes;
1022 }
1023 
1025 {
1026  return mNewZValues;
1027 }
1028 
1030 {
1031  return mNewXYValues;
1032 }
1033 
1035 {
1036  return mOldXYValues;
1037 }
1038 
1040 {
1041  return mNativeFacesIndexesGeometryChanged;
1042 }
1043 
1045 {
1046  return ( mFaceIndexesToRemove.isEmpty() &&
1047  mFacesToAdd.isEmpty() &&
1048  mFacesNeighborhoodToAdd.isEmpty() &&
1049  mFacesToRemove.isEmpty() &&
1050  mFacesNeighborhoodToRemove.isEmpty() &&
1051  mNeighborhoodChanges.isEmpty() &&
1052  mVerticesToAdd.isEmpty() &&
1053  mVertexToFaceToAdd.isEmpty() &&
1054  mVerticesToRemoveIndexes.isEmpty() &&
1055  mRemovedVertices.isEmpty() &&
1056  mVerticesToFaceRemoved.isEmpty() &&
1057  mVerticesToFaceChanges.isEmpty() &&
1058  mChangeCoordinateVerticesIndexes.isEmpty() &&
1059  mNewZValues.isEmpty() &&
1060  mOldZValues.isEmpty() &&
1061  mNewXYValues.isEmpty() &&
1062  mOldXYValues.isEmpty() &&
1063  mNativeFacesIndexesGeometryChanged.isEmpty() );
1064 }
1065 
1067 {
1068  return mVerticesToRemoveIndexes;
1069 }
1070 
1071 int QgsTopologicalMesh::Changes::addedFaceIndexInMesh( int internalIndex ) const
1072 {
1073  if ( internalIndex == -1 )
1074  return -1;
1075 
1076  return internalIndex + mAddedFacesFirstIndex;
1077 }
1078 
1079 int QgsTopologicalMesh::Changes::removedFaceIndexInMesh( int internalIndex ) const
1080 {
1081  if ( internalIndex == -1 )
1082  return -1;
1083 
1084  return mFaceIndexesToRemove.at( internalIndex );
1085 }
1086 
1088 {
1089  mAddedFacesFirstIndex = 0;
1090  mFaceIndexesToRemove.clear();
1091  mFacesToAdd.clear();
1092  mFacesNeighborhoodToAdd.clear();
1093  mFacesToRemove.clear();
1094  mFacesNeighborhoodToRemove.clear();
1095  mNeighborhoodChanges.clear();
1096 
1097  mVerticesToAdd.clear();
1098  mVertexToFaceToAdd.clear();
1099  mVerticesToRemoveIndexes.clear();
1100  mRemovedVertices.clear();
1101  mVerticesToFaceRemoved.clear();
1102  mVerticesToFaceChanges.clear();
1103 
1104  mChangeCoordinateVerticesIndexes.clear();
1105  mNewZValues.clear();
1106  mOldZValues.clear();
1107  mNewXYValues.clear();
1108  mOldXYValues.clear();
1109  mNativeFacesIndexesGeometryChanged.clear();
1110 }
1111 
1113 {
1114  Changes changes;
1115  changes.mVerticesToAdd.append( vertex );
1116  changes.mVertexToFaceToAdd.append( -1 );
1117 
1118  mMesh->vertices.append( vertex );
1119  mVertexToFace.append( -1 );
1120  referenceAsFreeVertex( mMesh->vertices.count() - 1 );
1121 
1122  return changes;
1123 }
1124 
1125 // Returns the orientation of the polygon formed by mesh vertices, <0 counter clockwise; >0 clockwise
1126 static double vertexPolygonOrientation( const QgsMesh &mesh, const QList<int> &vertexIndexes )
1127 {
1128  if ( vertexIndexes.count() < 3 )
1129  return 0.0;
1130 
1131  int hullDomainVertexPos = -1;
1132  double xMin = std::numeric_limits<double>::max();
1133  double yMin = std::numeric_limits<double>::max();
1134  for ( int i = 0; i < vertexIndexes.count(); ++i )
1135  {
1136  const QgsMeshVertex &vertex = mesh.vertices.at( vertexIndexes.at( i ) );
1137  if ( xMin >= vertex.x() && yMin > vertex.y() )
1138  {
1139  hullDomainVertexPos = i;
1140  xMin = vertex.x();
1141  yMin = vertex.y();
1142  }
1143  }
1144 
1145  if ( hullDomainVertexPos >= 0 )
1146  {
1147  int iv1 = vertexIndexes.at( ( hullDomainVertexPos - 1 + vertexIndexes.count() ) % vertexIndexes.count() );
1148  int iv2 = vertexIndexes.at( ( hullDomainVertexPos + 1 ) % vertexIndexes.count() );
1149  int ivc = vertexIndexes.at( ( hullDomainVertexPos ) );
1150  double cp = crossProduct( ivc, iv1, iv2, mesh );
1151  return cp;
1152  }
1153 
1154  return 0.0;
1155 }
1156 
1158 {
1159  if ( vertexIndex >= mVertexToFace.count() )
1160  return Changes();
1161 
1162  if ( mVertexToFace.at( vertexIndex ) == -1 ) //it is a free vertex
1163  {
1164  Changes changes;
1165  changes.mRemovedVertices.append( mMesh->vertices.at( vertexIndex ) );
1166  changes.mVerticesToRemoveIndexes.append( vertexIndex );
1167  changes.mVerticesToFaceRemoved.append( -1 );
1168  dereferenceAsFreeVertex( vertexIndex );
1169  mMesh->vertices[vertexIndex] = QgsMeshVertex();
1170  return changes;
1171  }
1172 
1173  //search concerned faces
1174  QgsMeshVertexCirculator circulator = vertexCirculator( vertexIndex );
1175  circulator.goBoundaryClockwise();
1176  QList<int> boundariesVertexIndex;
1177  QList<int> associateFaceToBoundaries;
1178  QList<int> removedFacesIndexes;
1179  QSet<int> boundaryInGlobalMesh;
1180 
1181  do
1182  {
1183  removedFacesIndexes.append( circulator.currentFaceIndex() );
1184  boundariesVertexIndex.append( circulator.oppositeVertexClockwise() );
1185  Q_ASSERT( !mMesh->vertices.at( boundariesVertexIndex.last() ).isEmpty() );
1186  const QgsMeshFace &currentFace = circulator.currentFace();
1187  associateFaceToBoundaries.append( mFacesNeighborhood.at( circulator.currentFaceIndex() ).at(
1188  vertexPositionInFace( boundariesVertexIndex.last(), currentFace ) ) );
1189 
1190  if ( currentFace.count() > 3 ) // quad or more, need other vertices
1191  {
1192  int posInface = vertexPositionInFace( vertexIndex, currentFace );
1193  for ( int i = 2; i < currentFace.count() - 1; ++i )
1194  {
1195  boundariesVertexIndex.append( currentFace.at( ( posInface + i ) % currentFace.count() ) );
1196  Q_ASSERT( !mMesh->vertices.at( boundariesVertexIndex.last() ).isEmpty() );
1197  associateFaceToBoundaries.append( mFacesNeighborhood.at( circulator.currentFaceIndex() ).at(
1198  vertexPositionInFace( boundariesVertexIndex.last(), currentFace ) ) );
1199  }
1200  }
1201  }
1202  while ( circulator.turnCounterClockwise() != -1 && circulator.currentFaceIndex() != removedFacesIndexes.first() );
1203 
1204  bool boundaryFill = false;
1205  if ( circulator.currentFaceIndex() == -1 ) //we are on boundary of the mesh
1206  {
1207  boundaryFill = true;
1208  //we need to add last vertex/boundary faces that was not added because we are on mesh boundary
1209  circulator.goBoundaryCounterClockwise();
1210  int lastVertexIndex = circulator.oppositeVertexCounterClockwise();
1211  boundariesVertexIndex.append( lastVertexIndex );
1212 
1213  // but we can be on the case where the last vertex share an edge with the first point,
1214  // in the case, the associate face on boundarie is not -1, but the face on the other side of the edge
1215  QgsMeshVertexCirculator boundaryCirculator = vertexCirculator( lastVertexIndex );
1216  boundaryCirculator.goBoundaryCounterClockwise();
1217  if ( boundaryCirculator.oppositeVertexCounterClockwise() == boundariesVertexIndex.first() )
1218  {
1219  associateFaceToBoundaries.append( boundaryCirculator.currentFaceIndex() );
1220  boundaryFill = false; //we are not a boundary fill anymore
1221  }
1222  else
1223  associateFaceToBoundaries.append( -1 );
1224 
1225  for ( const int index : std::as_const( boundariesVertexIndex ) )
1226  {
1227  if ( isVertexOnBoundary( index ) )
1228  boundaryInGlobalMesh.insert( index );
1229  }
1230  }
1231 
1232  int currentVertexToFace = mVertexToFace.at( vertexIndex );
1233  // here, we use the method removeFaces that effectivly removes and then constructs changes
1234  Changes changes = removeFaces( removedFacesIndexes );
1235 
1236  QList<QList<int>> holes;
1237  QList<QList<int>> associateMeshFacesToHoles;
1238 
1239  bool cancelOperation = false;
1240 
1241  if ( boundaryFill )
1242  {
1243  // The hole is not a closed polygon, we need to close it, but the closing segment can intersect another segments/vertices.
1244  // In this case we consider as many polygons as necessary.
1245 
1246  int startPos = 0;
1247  int finalPos = boundariesVertexIndex.count() - 1;
1248  QList<int> uncoveredVertex;
1249 
1250  QList<int> partToCheck = boundariesVertexIndex.mid( startPos, finalPos - startPos + 1 );
1251  QList<int> associateFacePart = associateFaceToBoundaries.mid( startPos, finalPos - startPos + 1 );
1252  while ( startPos < finalPos && !partToCheck.isEmpty() )
1253  {
1254  // check if we intersect an edge between first and second
1255  int secondPos = partToCheck.count() - 1;
1256  const QgsPoint &closingSegmentExtremety1 = mMesh->vertex( partToCheck.at( 0 ) );
1257  const QgsPoint &closingSegmentExtremety2 = mMesh->vertex( partToCheck.last() );
1258  bool isEdgeIntersect = false;
1259  for ( int i = 1; i < secondPos - 1; ++i )
1260  {
1261  const QgsPoint &p1 = mMesh->vertex( partToCheck.at( i ) );
1262  const QgsPoint &p2 = mMesh->vertex( partToCheck.at( i + 1 ) );
1263  bool isLineIntersection;
1264  QgsPoint intersectPoint;
1265  isEdgeIntersect = QgsGeometryUtils::segmentIntersection( closingSegmentExtremety1, closingSegmentExtremety2, p1, p2, intersectPoint, isLineIntersection, 1e-8, true );
1266  if ( isEdgeIntersect )
1267  break;
1268  }
1269 
1270  int index = partToCheck.at( 0 );
1271  if ( boundaryInGlobalMesh.contains( index ) && index != boundariesVertexIndex.at( 0 ) )
1272  {
1273  cancelOperation = true;
1274  break;
1275  }
1276 
1277  // if uncovered vertex is a boundary vertex in the global mesh (except first that is always a boundary in the global mesh)
1278  // This operation will leads to a unique shared vertex that is not allowed, you have to cancel the operation
1279  if ( isEdgeIntersect || vertexPolygonOrientation( *mMesh, partToCheck ) >= 0 )
1280  {
1281  partToCheck.removeLast();
1282  associateFacePart.removeAt( associateFacePart.count() - 2 );
1283  if ( partToCheck.count() == 1 )
1284  {
1285  uncoveredVertex.append( index );
1286  startPos = startPos + 1;
1287  partToCheck = boundariesVertexIndex.mid( startPos, finalPos - startPos + 1 );
1288  associateFacePart = associateFaceToBoundaries.mid( startPos, finalPos - startPos + 1 );
1289  }
1290  }
1291  else
1292  {
1293  // store the well defined hole
1294  holes.append( partToCheck );
1295  associateMeshFacesToHoles.append( associateFacePart );
1296 
1297  startPos = startPos + partToCheck.count() - 1;
1298  uncoveredVertex.append( partToCheck.at( 0 ) );
1299  partToCheck = boundariesVertexIndex.mid( startPos, finalPos - startPos + 1 );
1300  associateFacePart = associateFaceToBoundaries.mid( startPos, finalPos - startPos + 1 );
1301  }
1302  }
1303  }
1304  else
1305  {
1306  holes.append( boundariesVertexIndex );
1307  associateMeshFacesToHoles.append( associateFaceToBoundaries );
1308  }
1309 
1310  if ( cancelOperation )
1311  {
1312  reverseChanges( changes );
1313  return Changes();
1314  }
1315 
1316  Q_ASSERT( holes.count() == associateMeshFacesToHoles.count() );
1317 
1318  changes.mRemovedVertices.append( mMesh->vertices.at( vertexIndex ) );
1319  changes.mVerticesToRemoveIndexes.append( vertexIndex );
1320  changes.mVerticesToFaceRemoved.append( currentVertexToFace );
1321  // these changes contain information that will lead to reference the removed vertex as free vertex when reverse/reapply
1322  dereferenceAsFreeVertex( vertexIndex );
1323  mMesh->vertices[vertexIndex] = QgsMeshVertex();
1324  mVertexToFace[vertexIndex] = -1;
1325 
1326  int oldFacesCount = mMesh->faceCount();
1327  for ( int h = 0; h < holes.count(); ++h )
1328  {
1329  const QList<int> &holeVertices = holes.at( h );
1330  const QList<int> &associateMeshFacesToHole = associateMeshFacesToHoles.at( h );
1331  QHash<p2t::Point *, int> mapPoly2TriPointToVertex;
1332  std::vector<p2t::Point *> holeToFill( holeVertices.count() );
1333  try
1334  {
1335  for ( int i = 0; i < holeVertices.count(); ++i )
1336  {
1337  const QgsMeshVertex &vertex = mMesh->vertex( holeVertices.at( i ) );
1338  holeToFill[i] = new p2t::Point( vertex.x(), vertex.y() );
1339  mapPoly2TriPointToVertex.insert( holeToFill[i], holeVertices.at( i ) );
1340  }
1341 
1342  std::unique_ptr<p2t::CDT> cdt( new p2t::CDT( holeToFill ) );
1343 
1344  cdt->Triangulate();
1345  std::vector<p2t::Triangle *> triangles = cdt->GetTriangles();
1346  QVector<QgsMeshFace> newFaces( triangles.size() );
1347  for ( size_t i = 0; i < triangles.size(); ++i )
1348  {
1349  QgsMeshFace &face = newFaces[i];
1350  face.resize( 3 );
1351  for ( int j = 0; j < 3; j++ )
1352  {
1353  int vertInd = mapPoly2TriPointToVertex.value( triangles.at( i )->GetPoint( j ), -1 );
1354  if ( vertInd == -1 )
1355  throw std::exception();
1356  Q_ASSERT( !mMesh->vertices.at( vertInd ).isEmpty() );
1357  face[j] = vertInd;
1358  }
1359  }
1360 
1361  QgsMeshEditingError error;
1362  QgsTopologicalMesh::TopologicalFaces topologicalFaces = createNewTopologicalFaces( newFaces, false, error );
1364  throw std::exception();
1365  int newFaceIndexStartIndex = mMesh->faceCount();
1366  QgsTopologicalMesh::Changes addChanges;
1367  addChanges.mAddedFacesFirstIndex = newFaceIndexStartIndex;
1368  addChanges.mFacesToAdd = topologicalFaces.meshFaces();
1369  addChanges.mFacesNeighborhoodToAdd = topologicalFaces.mFacesNeighborhood;
1370 
1371  // vertices to face changes
1372  const QList<int> &verticesToFaceToChange = topologicalFaces.mVerticesToFace.keys();
1373  for ( const int vtc : verticesToFaceToChange )
1374  if ( mVertexToFace.at( vtc ) == -1 )
1375  addChanges.mVerticesToFaceChanges.append( {vtc,
1376  mVertexToFace.at( vtc ),
1377  addChanges.addedFaceIndexInMesh( topologicalFaces.mVerticesToFace.values( vtc ).first() ) } );
1378 
1379 
1380  // reindex neighborhood for new faces
1381  for ( int i = 0; i < topologicalFaces.mFaces.count(); ++i )
1382  {
1383  FaceNeighbors &faceNeighbors = addChanges.mFacesNeighborhoodToAdd[i];
1384  faceNeighbors = topologicalFaces.mFacesNeighborhood.at( i );
1385  for ( int n = 0; n < faceNeighbors.count(); ++n )
1386  {
1387  if ( faceNeighbors.at( n ) != -1 )
1388  faceNeighbors[n] += newFaceIndexStartIndex; //reindex internal neighborhood
1389  }
1390  }
1391 
1392  // link neighborhood for boundaries of each side
1393  for ( int i = 0 ; i < holeVertices.count(); ++i )
1394  {
1395  int vertexHoleIndex = holeVertices.at( i );
1396  int meshFaceBoundaryIndex = associateMeshFacesToHole.at( i );
1397 
1398  const QgsMeshVertexCirculator circulator = QgsMeshVertexCirculator( topologicalFaces, vertexHoleIndex );
1399  circulator.goBoundaryClockwise();
1400  int newFaceBoundaryLocalIndex = circulator.currentFaceIndex();
1401  int newFaceBoundaryIndexInMesh = circulator.currentFaceIndex() + newFaceIndexStartIndex;
1402  const QgsMeshFace &newFace = circulator.currentFace();
1403  int positionInNewFaces = vertexPositionInFace( vertexHoleIndex, newFace );
1404 
1405  if ( meshFaceBoundaryIndex != -1 )
1406  {
1407  const QgsMeshFace meshFace = mMesh->face( meshFaceBoundaryIndex );
1408  int positionInMeshFaceBoundary = vertexPositionInFace( *mMesh, vertexHoleIndex, meshFaceBoundaryIndex );
1409  positionInMeshFaceBoundary = ( positionInMeshFaceBoundary - 1 + meshFace.count() ) % meshFace.count(); //take the position just before
1410 
1411  addChanges.mNeighborhoodChanges.append( {meshFaceBoundaryIndex, positionInMeshFaceBoundary, -1, newFaceBoundaryIndexInMesh} );
1412  }
1413 
1414  addChanges.mFacesNeighborhoodToAdd[newFaceBoundaryLocalIndex][positionInNewFaces] = meshFaceBoundaryIndex;
1415  }
1416 
1417  applyChanges( addChanges );
1418 
1419  changes.mFacesToAdd.append( addChanges.mFacesToAdd );
1420  changes.mFacesNeighborhoodToAdd.append( addChanges.mFacesNeighborhoodToAdd );
1421  //for each neighborhood change, check if a corresponding change already exist and merge them, if not just append
1422  for ( const std::array<int, 4> &neighborChangeToAdd : std::as_const( addChanges.mNeighborhoodChanges ) )
1423  {
1424  bool merged = false;
1425  for ( std::array<int, 4> &existingNeighborChange : changes.mNeighborhoodChanges )
1426  {
1427  if ( existingNeighborChange.at( 0 ) == neighborChangeToAdd.at( 0 ) &&
1428  existingNeighborChange.at( 1 ) == neighborChangeToAdd.at( 1 ) )
1429  {
1430  merged = true;
1431  Q_ASSERT( existingNeighborChange.at( 3 ) == neighborChangeToAdd.at( 2 ) );
1432  existingNeighborChange[3] = neighborChangeToAdd.at( 3 );
1433  }
1434  }
1435  if ( !merged )
1436  changes.mNeighborhoodChanges.append( neighborChangeToAdd );
1437  }
1438  //for each vertex to face change, check if a corresponding change already exist and merge them, if not just append
1439  for ( const std::array<int, 3> &verticesToFaceToAdd : std::as_const( addChanges.mVerticesToFaceChanges ) )
1440  {
1441  bool merged = false;
1442  for ( std::array<int, 3> &existingVerticesToFace : changes.mVerticesToFaceChanges )
1443  {
1444  if ( existingVerticesToFace.at( 0 ) == verticesToFaceToAdd.at( 0 ) )
1445  {
1446  merged = true;
1447  Q_ASSERT( existingVerticesToFace.at( 2 ) == verticesToFaceToAdd.at( 1 ) );
1448  existingVerticesToFace[2] = verticesToFaceToAdd.at( 2 );
1449  }
1450  }
1451  if ( !merged )
1452  changes.mVerticesToFaceChanges.append( verticesToFaceToAdd );
1453  }
1454 
1455  qDeleteAll( holeToFill );
1456  }
1457  catch ( ... )
1458  {
1459  qDeleteAll( holeToFill );
1460  QgsMessageLog::logMessage( QObject::tr( "Triangulation failed. Skipping holeā€¦" ), QObject::tr( "Mesh Editing" ) );
1461  }
1462  }
1463  changes.mAddedFacesFirstIndex = oldFacesCount;
1464 
1465 
1466 
1467  changes.mAddedFacesFirstIndex = oldFacesCount;
1468 
1469  return changes;
1470 }
1471 
1473 {
1474  QSet<int> facesIndex;
1475  //Search for associated faces
1476  for ( int vertexIndex : vertices )
1477  facesIndex.unite( qgis::listToSet( facesAroundVertex( vertexIndex ) ) );
1478 
1479  // remove the faces
1480  Changes changes = removeFaces( facesIndex.values() );
1481 
1482  // removes the vertices
1483  for ( int vertexIndex : vertices )
1484  {
1485  int currentVertexToFace = mVertexToFace.at( vertexIndex );
1486  // here, we use the method removeFaces that effectivly removes and then constructs changes
1487  changes.mRemovedVertices.append( mMesh->vertices.at( vertexIndex ) );
1488  changes.mVerticesToRemoveIndexes.append( vertexIndex );
1489  changes.mVerticesToFaceRemoved.append( currentVertexToFace );
1490  // these changes contain information that will lead to reference the removed vertex as free vertex when reverse/reapply
1491  dereferenceAsFreeVertex( vertexIndex );
1492  mMesh->vertices[vertexIndex] = QgsMeshVertex();
1493  mVertexToFace[vertexIndex] = -1;
1494  }
1495 
1496  return changes;
1497 }
1498 
1500 {
1501  QList<int> boundariesToCheckClockwiseInNewFaces = topologicFaces.mBoundaries;
1502  QList<std::array<int, 2>> boundariesToCheckCounterClockwiseInNewFaces; //couple boundary / associate face in topologicFaces.mVerticesToFace
1503  QList<int> uniqueSharedVertexBoundary;
1504 
1505 
1506  // Go through the boundary and search for opposite boundary vertex clockwise in new faces and compare
1507  // with boundary opposite vertices on the mesh
1508  // If, in the mesh, the opposite vertex counter clockwise is not the same , another check will be done
1509  // later with counter clockwise in new faces
1510  // If, in the mesh, the opposite vertex clockwise is the same, this is an error
1511  while ( !boundariesToCheckClockwiseInNewFaces.isEmpty() )
1512  {
1513  int boundary = boundariesToCheckClockwiseInNewFaces.takeLast();
1514 
1515  const QList<int> &linkedFaces = topologicFaces.mVerticesToFace.values( boundary );
1516 
1517  for ( int const linkedFace : linkedFaces )
1518  {
1519 
1520  //if the boundary is a free vertex in the destination mesh, no need to check
1521  if ( mVertexToFace.at( boundary ) == -1 )
1522  continue;
1523 
1524  QgsMeshVertexCirculator newFacescirculator( topologicFaces, linkedFace, boundary );
1525  QgsMeshVertexCirculator meshCirculator = vertexCirculator( boundary );
1526 
1527  if ( !newFacescirculator.isValid() )
1529 
1530  //Get the opposite vertex on the clockwise side with new faces block
1531  if ( !newFacescirculator.goBoundaryClockwise() )
1533 
1534  int oppositeVertexInNewFaces = newFacescirculator.oppositeVertexClockwise();
1535 
1536  if ( !meshCirculator.goBoundaryCounterClockwise() )
1538 
1539  int oppositeVertexCCWInMesh = meshCirculator.oppositeVertexCounterClockwise();
1540 
1541  if ( oppositeVertexCCWInMesh == oppositeVertexInNewFaces ) //this boundary is OK, continue wit next one
1542  continue;
1543  else
1544  {
1545  //to avoid manifold face that could pass through the check, compare not only the boundary edges but also with the opposite internal edge in new face
1546  const QgsMeshFace &newFaceOnBoundary = newFacescirculator.currentFace();
1547  int faceSize = newFaceOnBoundary.size();
1548  int posInNewFace = vertexPositionInFace( boundary, newFaceOnBoundary );
1549  int previousVertexIndex = ( posInNewFace + faceSize - 1 ) % faceSize;
1550  if ( newFaceOnBoundary.at( previousVertexIndex ) == oppositeVertexCCWInMesh )
1552  }
1553 
1554  meshCirculator.goBoundaryClockwise();
1555 
1556  int oppositeVertexCWInMesh = meshCirculator.oppositeVertexClockwise();
1557 
1558  if ( oppositeVertexCWInMesh == oppositeVertexInNewFaces )
1560 
1561  //if we are here we need more checks
1562  boundariesToCheckCounterClockwiseInNewFaces.append( {boundary, linkedFace} );
1563  }
1564  }
1565 
1566  // Check now with opposite boundary vertex counter clockwise in new faces
1567  while ( !boundariesToCheckCounterClockwiseInNewFaces.isEmpty() )
1568  {
1569  std::array<int, 2> boundaryLinkedface = boundariesToCheckCounterClockwiseInNewFaces.takeLast();
1570  int boundary = boundaryLinkedface.at( 0 );
1571  int linkedFace = boundaryLinkedface.at( 1 );
1572 
1573  QgsMeshVertexCirculator newFacescirculator( topologicFaces, linkedFace, boundary );
1574  QgsMeshVertexCirculator meshCirculator = vertexCirculator( boundary );
1575 
1576  if ( !newFacescirculator.goBoundaryCounterClockwise() )
1578 
1579  int oppositeVertexInNewFaces = newFacescirculator.oppositeVertexCounterClockwise();
1580 
1581  if ( !meshCirculator.goBoundaryClockwise() )
1583 
1584  int oppositeVertexCWInMesh = meshCirculator.oppositeVertexClockwise();
1585 
1586  if ( oppositeVertexCWInMesh == oppositeVertexInNewFaces ) //this boundary is OK, continue with next one
1587  continue;
1588 
1589  meshCirculator.goBoundaryCounterClockwise();
1590 
1591  int oppositeVertexCCWInMesh = meshCirculator.oppositeVertexCounterClockwise();
1592 
1593  if ( oppositeVertexCCWInMesh == oppositeVertexInNewFaces )
1595 
1596  uniqueSharedVertexBoundary.append( boundary );
1597  }
1598 
1599  if ( !uniqueSharedVertexBoundary.isEmpty() )
1600  return QgsMeshEditingError( Qgis::MeshEditingErrorType::UniqueSharedVertex, uniqueSharedVertexBoundary.first() );
1601 
1602  // Check if internal vertices of new faces block are free in the mesh
1603  QSet<int> boundaryVertices = qgis::listToSet( topologicFaces.mBoundaries );
1604  for ( const QgsMeshFace &newFace : std::as_const( topologicFaces.mFaces ) )
1605  {
1606  for ( const int vertexIndex : newFace )
1607  {
1608  if ( boundaryVertices.contains( vertexIndex ) )
1609  continue;
1610  if ( mVertexToFace.at( vertexIndex ) != -1 )
1612  }
1613  }
1614 
1615  return QgsMeshEditingError();
1616 }
1617 
1619 {
1620  mFaces.clear();
1621  mFacesNeighborhood.clear();
1622  mVerticesToFace.clear();
1623  mBoundaries.clear();
1624 }
1625 
1626 QVector<QgsTopologicalMesh::FaceNeighbors> QgsTopologicalMesh::TopologicalFaces::facesNeighborhood() const
1627 {
1628  return mFacesNeighborhood;
1629 }
1630 
1632 {
1633  if ( mVerticesToFace.contains( vertexIndex ) )
1634  return mVerticesToFace.values( vertexIndex ).at( 0 );
1635 
1636  return -1;
1637 }
1638 
1640 {
1641  QgsTopologicalMesh topologicMesh;
1642  topologicMesh.mMesh = mesh;
1643  topologicMesh.mVertexToFace = QVector<int>( mesh->vertexCount(), -1 );
1644  topologicMesh.mMaximumVerticesPerFace = maxVerticesPerFace;
1646 
1647  for ( int i = 0; i < mesh->faceCount(); ++i )
1648  {
1649  if ( mesh->face( i ).isEmpty() )
1650  continue;
1651  if ( maxVerticesPerFace != 0 && mesh->face( i ).count() > maxVerticesPerFace )
1652  {
1654  break;
1655  }
1656 
1657  error = counterClockwiseFaces( mesh->faces[i], mesh );
1659  {
1661  error.elementIndex = i;
1662  break;
1663  }
1664  }
1665 
1667  {
1668  TopologicalFaces subMesh = createTopologicalFaces( mesh->faces, &topologicMesh.mVertexToFace, error, false );
1669  topologicMesh.mFacesNeighborhood = subMesh.mFacesNeighborhood;
1670 
1671  for ( int i = 0; i < topologicMesh.mMesh->vertexCount(); ++i )
1672  {
1673  if ( topologicMesh.mVertexToFace.at( i ) == -1 )
1674  topologicMesh.mFreeVertices.insert( i );
1675  }
1676  }
1677 
1678  return topologicMesh;
1679 }
1680 
1681 QgsTopologicalMesh::TopologicalFaces QgsTopologicalMesh::createNewTopologicalFaces( const QVector<QgsMeshFace> &faces, bool uniqueSharedVertexAllowed, QgsMeshEditingError &error )
1682 {
1683  return createTopologicalFaces( faces, nullptr, error, uniqueSharedVertexAllowed );
1684 }
1685 
1686 
1687 QgsTopologicalMesh::TopologicalFaces QgsTopologicalMesh::createTopologicalFaces(
1688  const QVector<QgsMeshFace> &faces,
1689  QVector<int> *globalVertexToFace,
1690  QgsMeshEditingError &error,
1691  bool allowUniqueSharedVertex )
1692 {
1693  int facesCount = faces.count();
1694  QVector<FaceNeighbors> faceTopologies;
1695  QMultiHash<int, int> verticesToFace;
1696 
1697  error = QgsMeshEditingError();
1698  TopologicalFaces ret;
1699 
1700  // Contains for each vertex a map (opposite vertex # edge) --> face index
1701  // when turning counter clockwise, if v1 first vertex and v2 second one, [v1][v2]--> neighbor face
1702  QMap<int, QMap<int, int>> verticesToNeighbor;
1703 
1704  for ( int faceIndex = 0; faceIndex < facesCount; ++faceIndex )
1705  {
1706  const QgsMeshFace &face = faces.at( faceIndex );
1707  int faceSize = face.count();
1708  //Fill vertices to neighbor faces
1709  for ( int i = 0; i < faceSize; ++i )
1710  {
1711  int v1 = face[i % faceSize];
1712  int v2 = face[( i + 1 ) % faceSize];
1713  if ( verticesToNeighbor[v2].contains( v1 ) )
1714  {
1716  return ret;
1717  }
1718  else
1719  verticesToNeighbor[v2].insert( v1, faceIndex );
1720  }
1721  }
1722 
1723  faceTopologies = QVector<FaceNeighbors>( faces.count() );
1724 
1725  QSet<int> boundaryVertices;
1726 
1727  for ( int faceIndex = 0; faceIndex < facesCount; ++faceIndex )
1728  {
1729  const QgsMeshFace &face = faces.at( faceIndex );
1730  int faceSize = face.size();
1731  FaceNeighbors &faceTopology = faceTopologies[faceIndex];
1732  faceTopology.resize( faceSize );
1733 
1734  for ( int i = 0; i < faceSize; ++i )
1735  {
1736  int v1 = face.at( i );
1737  int v2 = face.at( ( i + 1 ) % faceSize );
1738 
1739  if ( globalVertexToFace )
1740  {
1741  if ( ( *globalVertexToFace )[v1] == -1 )
1742  ( *globalVertexToFace )[v1] = faceIndex ;
1743  }
1744  else
1745  {
1746  if ( allowUniqueSharedVertex || !verticesToFace.contains( v1 ) )
1747  verticesToFace.insert( v1, faceIndex ) ;
1748  }
1749 
1750  QMap<int, int> &edges = verticesToNeighbor[v1];
1751  if ( edges.contains( v2 ) )
1752  faceTopology[i] = edges.value( v2 );
1753  else
1754  {
1755  faceTopology[i] = -1;
1756 
1757  if ( !allowUniqueSharedVertex )
1758  {
1759  if ( boundaryVertices.contains( v1 ) )
1760  {
1761  error = QgsMeshEditingError( Qgis::MeshEditingErrorType::UniqueSharedVertex, v1 ); // if a vertices is more than one time in the boundary, that means faces share only one vertices
1762  return ret;
1763  }
1764  }
1765  boundaryVertices.insert( v1 );
1766  }
1767  }
1768  }
1769 
1770  ret.mFaces = faces;
1771  ret.mFacesNeighborhood = faceTopologies;
1772  ret.mBoundaries = boundaryVertices.values();
1773  ret.mVerticesToFace = verticesToFace;
1774  return ret;
1775 }
1776 
1777 QVector<int> QgsTopologicalMesh::neighborsOfFace( int faceIndex ) const
1778 {
1779  return mFacesNeighborhood.at( faceIndex );
1780 }
1781 
1782 QList<int> QgsTopologicalMesh::facesAroundVertex( int vertexIndex ) const
1783 {
1784  QgsMeshVertexCirculator circ = vertexCirculator( vertexIndex );
1785 
1786  return circ.facesAround();
1787 }
1788 
1790 {
1791  QSet<int> removedFaces = qgis::listToSet( facesIndexes );
1792  QSet<int> concernedFaces = concernedFacesBy( facesIndexes );
1793 
1794  for ( const int f : std::as_const( removedFaces ) )
1795  concernedFaces.remove( f );
1796 
1797  QVector<QgsMeshFace> remainingFaces;
1798  remainingFaces.reserve( concernedFaces.count() );
1799  for ( const int f : std::as_const( concernedFaces ) )
1800  remainingFaces.append( mMesh->face( f ) );
1801 
1802  QgsMeshEditingError error;
1803  createTopologicalFaces( remainingFaces, nullptr, error, false );
1804 
1805  return error;
1806 }
1807 
1808 QgsTopologicalMesh::Changes QgsTopologicalMesh::removeFaces( const QList<int> facesIndexesToRemove )
1809 {
1810  Changes changes;
1811  changes.mFaceIndexesToRemove = facesIndexesToRemove;
1812  changes.mFacesToRemove.resize( facesIndexesToRemove.count() );
1813  changes.mFacesNeighborhoodToRemove.resize( facesIndexesToRemove.count() );
1814 
1815  QSet<int> indexSet = qgis::listToSet( facesIndexesToRemove );
1816  QSet<int> threatedVertex;
1817 
1818  for ( int i = 0; i < facesIndexesToRemove.count(); ++i )
1819  {
1820  const int faceIndex = facesIndexesToRemove.at( i );
1821  const QgsMeshFace &face = mMesh->face( faceIndex );
1822  changes.mFacesToRemove[i] = face;
1823  const FaceNeighbors &neighborhood = mFacesNeighborhood.at( faceIndex );
1824  changes.mFacesNeighborhoodToRemove[i] = neighborhood;
1825  for ( int j = 0; j < face.count(); ++j )
1826  {
1827  //change the neighborhood for each neighbor face
1828  int neighborIndex = neighborhood.at( j );
1829  if ( neighborIndex != -1 && !indexSet.contains( neighborIndex ) )
1830  {
1831  int positionInNeighbor = mFacesNeighborhood.at( neighborIndex ).indexOf( faceIndex );
1832  changes.mNeighborhoodChanges.append( {neighborIndex, positionInNeighbor, faceIndex, -1} );
1833  }
1834 
1835  //change vertexToFace
1836  int vertexIndex = face.at( j );
1837  if ( !threatedVertex.contains( vertexIndex ) && indexSet.contains( mVertexToFace.at( vertexIndex ) ) )
1838  {
1839  int oldValue = mVertexToFace.at( vertexIndex );
1840  //look for another face linked to this vertex
1841  int refValue = -1;
1842  if ( neighborIndex != -1 && !indexSet.contains( neighborIndex ) ) //if exist, simpler to take it
1843  refValue = neighborIndex;
1844  else
1845  {
1846  QList<int> aroundFaces = facesAroundVertex( vertexIndex );
1847  aroundFaces.removeOne( faceIndex );
1848  if ( !aroundFaces.isEmpty() )
1849  {
1850  while ( !aroundFaces.isEmpty() && refValue == -1 )
1851  {
1852  if ( !indexSet.contains( aroundFaces.first() ) )
1853  refValue = aroundFaces.first();
1854  else
1855  aroundFaces.removeFirst();
1856  }
1857  }
1858  }
1859  changes.mVerticesToFaceChanges.append( {vertexIndex, oldValue, refValue} );
1860  threatedVertex.insert( vertexIndex );
1861  }
1862  }
1863  }
1864 
1865  applyChanges( changes );
1866 
1867  return changes;
1868 }
1869 
1870 bool QgsTopologicalMesh::eitherSideFacesAndVertices( int vertexIndex1,
1871  int vertexIndex2,
1872  int &face1,
1873  int &face2,
1874  int &neighborVertex1InFace1,
1875  int &neighborVertex1InFace2,
1876  int &neighborVertex2inFace1,
1877  int &neighborVertex2inFace2 ) const
1878 {
1879  QgsMeshVertexCirculator circulator1 = vertexCirculator( vertexIndex1 );
1880  QgsMeshVertexCirculator circulator2 = vertexCirculator( vertexIndex2 );
1881 
1882  circulator1.goBoundaryClockwise();
1883  int firstFace1 = circulator1.currentFaceIndex();
1884  circulator2.goBoundaryClockwise();
1885  int firstFace2 = circulator2.currentFaceIndex();
1886 
1887  if ( circulator1.oppositeVertexCounterClockwise() != vertexIndex2 )
1888  while ( circulator1.turnCounterClockwise() != -1 &&
1889  circulator1.currentFaceIndex() != firstFace1 &&
1890  circulator1.oppositeVertexCounterClockwise() != vertexIndex2 ) {}
1891 
1892  if ( circulator2.oppositeVertexCounterClockwise() != vertexIndex1 )
1893  while ( circulator2.turnCounterClockwise() != -1 &&
1894  circulator2.currentFaceIndex() != firstFace2 &&
1895  circulator2.oppositeVertexCounterClockwise() != vertexIndex1 ) {}
1896 
1897  if ( circulator1.oppositeVertexCounterClockwise() != vertexIndex2
1898  || circulator2.oppositeVertexCounterClockwise() != vertexIndex1 )
1899  return false;
1900 
1901  face1 = circulator1.currentFaceIndex();
1902  face2 = circulator2.currentFaceIndex();
1903 
1904  neighborVertex1InFace1 = circulator1.oppositeVertexClockwise();
1905  circulator1.turnCounterClockwise();
1906  neighborVertex1InFace2 = circulator1.oppositeVertexCounterClockwise();
1907 
1908  neighborVertex2inFace2 = circulator2.oppositeVertexClockwise();
1909  circulator2.turnCounterClockwise();
1910  neighborVertex2inFace1 = circulator2.oppositeVertexCounterClockwise();
1911 
1912  return true;
1913 }
1914 
1915 bool QgsTopologicalMesh::edgeCanBeFlipped( int vertexIndex1, int vertexIndex2 ) const
1916 {
1917  int faceIndex1;
1918  int faceIndex2;
1919  int oppositeVertexFace1;
1920  int oppositeVertexFace2;
1921  int supposedOppositeVertexFace1;
1922  int supposedoppositeVertexFace2;
1923 
1924  bool result = eitherSideFacesAndVertices(
1925  vertexIndex1,
1926  vertexIndex2,
1927  faceIndex1,
1928  faceIndex2,
1929  oppositeVertexFace1,
1930  supposedoppositeVertexFace2,
1931  supposedOppositeVertexFace1,
1932  oppositeVertexFace2 );
1933 
1934  if ( !result ||
1935  faceIndex1 < 0 ||
1936  faceIndex2 < 0 ||
1937  oppositeVertexFace1 < 0 ||
1938  oppositeVertexFace2 < 0 ||
1939  supposedOppositeVertexFace1 != oppositeVertexFace1 ||
1940  supposedoppositeVertexFace2 != oppositeVertexFace2 )
1941  return false;
1942 
1943  const QgsMeshFace &face1 = mMesh->face( faceIndex1 );
1944  const QgsMeshFace &face2 = mMesh->face( faceIndex2 );
1945 
1946 
1947  if ( face1.count() != 3 || face2.count() != 3 )
1948  return false;
1949 
1950  double crossProduct1 = crossProduct( vertexIndex1, oppositeVertexFace1, oppositeVertexFace2, *mMesh );
1951  double crossProduct2 = crossProduct( vertexIndex2, oppositeVertexFace1, oppositeVertexFace2, *mMesh );
1952 
1953  return crossProduct1 * crossProduct2 < 0;
1954 }
1955 
1956 QgsTopologicalMesh::Changes QgsTopologicalMesh::flipEdge( int vertexIndex1, int vertexIndex2 )
1957 {
1958  int faceIndex1;
1959  int faceIndex2;
1960  int oppositeVertexFace1;
1961  int oppositeVertexFace2;
1962  int supposedOppositeVertexFace1;
1963  int supposedoppositeVertexFace2;
1964 
1965  bool result = eitherSideFacesAndVertices(
1966  vertexIndex1,
1967  vertexIndex2,
1968  faceIndex1,
1969  faceIndex2,
1970  oppositeVertexFace1,
1971  supposedoppositeVertexFace2,
1972  supposedOppositeVertexFace1,
1973  oppositeVertexFace2 );
1974 
1975  if ( !result ||
1976  faceIndex1 < 0 ||
1977  faceIndex2 < 0 ||
1978  oppositeVertexFace1 < 0 ||
1979  oppositeVertexFace2 < 0 ||
1980  supposedOppositeVertexFace1 != oppositeVertexFace1 ||
1981  supposedoppositeVertexFace2 != oppositeVertexFace2 )
1982  return Changes();
1983 
1984 
1985  Changes changes;
1986 
1987  const QgsMeshFace &face1 = mMesh->face( faceIndex1 );
1988  const QgsMeshFace &face2 = mMesh->face( faceIndex2 );
1989 
1990  Q_ASSERT( face1.count() == 3 );
1991  Q_ASSERT( face2.count() == 3 );
1992 
1993  int pos1 = vertexPositionInFace( vertexIndex1, face1 );
1994  int pos2 = vertexPositionInFace( vertexIndex2, face2 );
1995 
1996  int neighborFace1 = mFacesNeighborhood.at( faceIndex1 ).at( pos1 );
1997  int posInNeighbor1 = vertexPositionInFace( *mMesh, oppositeVertexFace1, neighborFace1 );
1998  int neighborFace2 = mFacesNeighborhood.at( faceIndex1 ).at( ( pos1 + 1 ) % 3 );
1999  int posInNeighbor2 = vertexPositionInFace( *mMesh, vertexIndex2, neighborFace2 );
2000  int neighborFace3 = mFacesNeighborhood.at( faceIndex2 ).at( pos2 );
2001  int posInNeighbor3 = vertexPositionInFace( *mMesh, oppositeVertexFace2, neighborFace3 );
2002  int neighborFace4 = mFacesNeighborhood.at( faceIndex2 ).at( ( pos2 + 1 ) % 3 );
2003  int posInNeighbor4 = vertexPositionInFace( *mMesh, vertexIndex1, neighborFace4 );
2004 
2005  changes.mFaceIndexesToRemove.append( faceIndex1 );
2006  changes.mFaceIndexesToRemove.append( faceIndex2 );
2007  changes.mFacesToRemove.append( face1 );
2008  changes.mFacesToRemove.append( face2 );
2009  changes.mFacesNeighborhoodToRemove.append( mFacesNeighborhood.at( faceIndex1 ) );
2010  changes.mFacesNeighborhoodToRemove.append( mFacesNeighborhood.at( faceIndex2 ) );
2011  int startIndex = mMesh->faceCount();
2012  changes.mAddedFacesFirstIndex = startIndex;
2013  changes.mFacesToAdd.append( {oppositeVertexFace1, oppositeVertexFace2, vertexIndex1} );
2014  changes.mFacesToAdd.append( {oppositeVertexFace2, oppositeVertexFace1, vertexIndex2} );
2015  changes.mFacesNeighborhoodToAdd.append( {startIndex + 1,
2016  mFacesNeighborhood.at( faceIndex2 ).at( ( pos2 + 1 ) % 3 ),
2017  mFacesNeighborhood.at( faceIndex1 ).at( pos1 )} );
2018  changes.mFacesNeighborhoodToAdd.append( {startIndex,
2019  mFacesNeighborhood.at( faceIndex1 ).at( ( pos1 + 1 ) % 3 ),
2020  mFacesNeighborhood.at( faceIndex2 ).at( pos2 )} );
2021 
2022  if ( neighborFace1 >= 0 )
2023  changes.mNeighborhoodChanges.append( {neighborFace1, posInNeighbor1, faceIndex1, startIndex} );
2024  if ( neighborFace2 >= 0 )
2025  changes.mNeighborhoodChanges.append( {neighborFace2, posInNeighbor2, faceIndex1, startIndex + 1} );
2026  if ( neighborFace3 >= 0 )
2027  changes.mNeighborhoodChanges.append( {neighborFace3, posInNeighbor3, faceIndex2, startIndex + 1} );
2028  if ( neighborFace4 >= 0 )
2029  changes.mNeighborhoodChanges.append( {neighborFace4, posInNeighbor4, faceIndex2, startIndex} );
2030 
2031 
2032  if ( mVertexToFace.at( vertexIndex1 ) == faceIndex1 || mVertexToFace.at( vertexIndex1 ) == faceIndex2 )
2033  changes.mVerticesToFaceChanges.append( {vertexIndex1, mVertexToFace.at( vertexIndex1 ), startIndex} );
2034  if ( mVertexToFace.at( vertexIndex2 ) == faceIndex1 || mVertexToFace.at( vertexIndex2 ) == faceIndex2 )
2035  changes.mVerticesToFaceChanges.append( {vertexIndex2, mVertexToFace.at( vertexIndex2 ), startIndex + 1} );
2036 
2037  if ( mVertexToFace.at( oppositeVertexFace1 ) == faceIndex1 )
2038  changes.mVerticesToFaceChanges.append( {oppositeVertexFace1, faceIndex1, startIndex} );
2039 
2040  if ( mVertexToFace.at( oppositeVertexFace2 ) == faceIndex2 )
2041  changes.mVerticesToFaceChanges.append( {oppositeVertexFace2, faceIndex2, startIndex + 1} );
2042 
2043  applyChanges( changes );
2044 
2045  return changes;
2046 }
2047 
2048 bool QgsTopologicalMesh::canBeMerged( int vertexIndex1, int vertexIndex2 ) const
2049 {
2050  int faceIndex1;
2051  int faceIndex2;
2052  int neighborVertex1InFace1;
2053  int neighborVertex1InFace2;
2054  int neighborVertex2inFace1;
2055  int neighborVertex2inFace2;
2056 
2057  bool result = eitherSideFacesAndVertices(
2058  vertexIndex1,
2059  vertexIndex2,
2060  faceIndex1,
2061  faceIndex2,
2062  neighborVertex1InFace1,
2063  neighborVertex1InFace2,
2064  neighborVertex2inFace1,
2065  neighborVertex2inFace2 );
2066 
2067  if ( !result ||
2068  faceIndex1 < 0 ||
2069  faceIndex2 < 0 )
2070  return false;
2071 
2072  const QgsMeshFace &face1 = mMesh->face( faceIndex1 );
2073  const QgsMeshFace &face2 = mMesh->face( faceIndex2 );
2074 
2075  if ( face1.count() + face2.count() - 2 > mMaximumVerticesPerFace )
2076  return false;
2077 
2078  QgsMeshVertex v1 = mMesh->vertices.at( vertexIndex1 );
2079  QgsMeshVertex v2 = mMesh->vertices.at( vertexIndex2 );
2080  QgsMeshVertex nv11 = mMesh->vertices.at( neighborVertex1InFace1 );
2081  QgsMeshVertex nv12 = mMesh->vertices.at( neighborVertex1InFace2 );
2082  QgsMeshVertex nv21 = mMesh->vertices.at( neighborVertex2inFace1 );
2083  QgsMeshVertex nv22 = mMesh->vertices.at( neighborVertex2inFace2 );
2084 
2085  double crossProduct1 = crossProduct( vertexIndex1, neighborVertex1InFace1, neighborVertex1InFace2, *mMesh );
2086  double crossProduct2 = crossProduct( vertexIndex2, neighborVertex2inFace1, neighborVertex2inFace2, *mMesh );
2087 
2088  return crossProduct1 * crossProduct2 < 0;
2089 }
2090 
2091 QgsTopologicalMesh::Changes QgsTopologicalMesh::merge( int vertexIndex1, int vertexIndex2 )
2092 {
2093  int faceIndex1;
2094  int faceIndex2;
2095  int neighborVertex1InFace1;
2096  int neighborVertex1InFace2;
2097  int neighborVertex2inFace1;
2098  int neighborVertex2inFace2;
2099 
2100  bool result = eitherSideFacesAndVertices(
2101  vertexIndex1,
2102  vertexIndex2,
2103  faceIndex1,
2104  faceIndex2,
2105  neighborVertex1InFace1,
2106  neighborVertex1InFace2,
2107  neighborVertex2inFace1,
2108  neighborVertex2inFace2 );
2109 
2110  if ( !result ||
2111  faceIndex1 < 0 ||
2112  faceIndex2 < 0 )
2113  return Changes();
2114 
2115  Changes changes;
2116 
2117  const QgsMeshFace &face1 = mMesh->face( faceIndex1 );
2118  const QgsMeshFace &face2 = mMesh->face( faceIndex2 );
2119  int faceSize1 = face1.count();
2120  int faceSize2 = face2.count();
2121 
2122  int pos1 = vertexPositionInFace( vertexIndex1, face1 );
2123  int pos2 = vertexPositionInFace( vertexIndex2, face2 );
2124 
2125  changes.mFaceIndexesToRemove.append( faceIndex1 );
2126  changes.mFaceIndexesToRemove.append( faceIndex2 );
2127  changes.mFacesToRemove.append( face1 );
2128  changes.mFacesToRemove.append( face2 );
2129  changes.mFacesNeighborhoodToRemove.append( mFacesNeighborhood.at( faceIndex1 ) );
2130  changes.mFacesNeighborhoodToRemove.append( mFacesNeighborhood.at( faceIndex2 ) );
2131  int startIndex = mMesh->faceCount();
2132  changes.mAddedFacesFirstIndex = startIndex;
2133 
2134  QgsMeshFace newface;
2135  FaceNeighbors newNeighborhood;
2136  for ( int i = 0; i < faceSize1 - 1; ++i )
2137  {
2138  int currentPos = ( pos1 + i ) % faceSize1;
2139  newface.append( face1.at( currentPos ) ); //add vertex of face1
2140 
2141  int currentNeighbor = mFacesNeighborhood.at( faceIndex1 ).at( currentPos );
2142  newNeighborhood.append( currentNeighbor );
2143 
2144  if ( currentNeighbor != -1 )
2145  {
2146  int currentPosInNeighbor = vertexPositionInFace( *mMesh, face1.at( ( currentPos + 1 ) % faceSize1 ), currentNeighbor );
2147  changes.mNeighborhoodChanges.append( {currentNeighbor, currentPosInNeighbor, faceIndex1, startIndex} );
2148  }
2149  }
2150  for ( int i = 0; i < faceSize2 - 1; ++i )
2151  {
2152  int currentPos = ( pos2 + i ) % faceSize2;
2153  newface.append( face2.at( currentPos ) ); //add vertex of face2
2154 
2155  int currentNeighbor = mFacesNeighborhood.at( faceIndex2 ).at( currentPos );
2156  newNeighborhood.append( currentNeighbor );
2157 
2158  if ( currentNeighbor != -1 )
2159  {
2160  int currentPosInNeighbor = vertexPositionInFace( *mMesh, face2.at( ( currentPos + 1 ) % faceSize2 ), currentNeighbor );
2161  changes.mNeighborhoodChanges.append( {currentNeighbor, currentPosInNeighbor, faceIndex2, startIndex} );
2162  }
2163  }
2164 
2165  for ( int i = 0; i < faceSize1; ++i )
2166  if ( mVertexToFace.at( face1.at( i ) ) == faceIndex1 )
2167  changes.mVerticesToFaceChanges.append( {face1.at( i ), faceIndex1, startIndex} );
2168 
2169  for ( int i = 0; i < faceSize2; ++i )
2170  if ( mVertexToFace.at( face2.at( i ) ) == faceIndex2 )
2171  changes.mVerticesToFaceChanges.append( {face2.at( i ), faceIndex2, startIndex} );
2172 
2173  changes.mFacesToAdd.append( newface );
2174  changes.mFacesNeighborhoodToAdd.append( newNeighborhood );
2175 
2176  applyChanges( changes );
2177 
2178  return changes;
2179 }
2180 
2181 bool QgsTopologicalMesh::canBeSplit( int faceIndex ) const
2182 {
2183  const QgsMeshFace face = mMesh->face( faceIndex );
2184 
2185  return face.count() == 4;
2186 }
2187 
2189 {
2190  //search for the spliited angle (greater angle)
2191  const QgsMeshFace &face = mMesh->face( faceIndex );
2192  int faceSize = face.count();
2193 
2194  Q_ASSERT( faceSize == 4 );
2195 
2196  double maxAngle = 0;
2197  int splitVertexPos = -1;
2198  for ( int i = 0; i < faceSize; ++i )
2199  {
2200  QgsVector vect1( mMesh->vertex( face.at( i ) ) - mMesh->vertex( face.at( ( i + 1 ) % faceSize ) ) );
2201  QgsVector vect2( mMesh->vertex( face.at( ( i + 2 ) % faceSize ) ) - mMesh->vertex( face.at( ( i + 1 ) % faceSize ) ) );
2202 
2203  double angle = std::abs( vect1.angle( vect2 ) );
2204  angle = std::min( angle, 2.0 * M_PI - angle );
2205  if ( angle > maxAngle )
2206  {
2207  maxAngle = angle;
2208  splitVertexPos = ( i + 1 ) % faceSize;
2209  }
2210  }
2211 
2212  Changes changes;
2213 
2214  const QgsMeshFace newFace1 = {face.at( splitVertexPos ),
2215  face.at( ( splitVertexPos + 1 ) % faceSize ),
2216  face.at( ( splitVertexPos + 2 ) % faceSize )
2217  };
2218 
2219  const QgsMeshFace newFace2 = {face.at( splitVertexPos ),
2220  face.at( ( splitVertexPos + 2 ) % faceSize ),
2221  face.at( ( splitVertexPos + 3 ) % faceSize )
2222  };
2223 
2224  QVector<int> neighborIndex( faceSize );
2225  QVector<int> posInNeighbor( faceSize );
2226 
2227  for ( int i = 0; i < faceSize; ++i )
2228  {
2229  neighborIndex[i] = mFacesNeighborhood.at( faceIndex ).at( ( splitVertexPos + i ) % faceSize );
2230  posInNeighbor[i] = vertexPositionInFace( *mMesh, face.at( ( splitVertexPos + i + 1 ) % faceSize ), neighborIndex[i] );
2231  }
2232 
2233  changes.mFaceIndexesToRemove.append( faceIndex );
2234  changes.mFacesToRemove.append( face );
2235  changes.mFacesNeighborhoodToRemove.append( mFacesNeighborhood.at( faceIndex ) );
2236  int startIndex = mMesh->faceCount();
2237  changes.mAddedFacesFirstIndex = startIndex;
2238  changes.mFacesToAdd.append( newFace1 );
2239  changes.mFacesToAdd.append( newFace2 );
2240 
2241  changes.mFacesNeighborhoodToAdd.append( {mFacesNeighborhood.at( faceIndex ).at( splitVertexPos ),
2242  mFacesNeighborhood.at( faceIndex ).at( ( splitVertexPos + 1 ) % faceSize ),
2243  startIndex + 1} );
2244  changes.mFacesNeighborhoodToAdd.append( {startIndex,
2245  mFacesNeighborhood.at( faceIndex ).at( ( splitVertexPos + 2 ) % faceSize ),
2246  mFacesNeighborhood.at( faceIndex ).at( ( splitVertexPos + 3 ) % faceSize )} );
2247 
2248  for ( int i = 0; i < faceSize; ++i )
2249  {
2250  if ( neighborIndex[i] >= 0 )
2251  changes.mNeighborhoodChanges.append( {neighborIndex[i], posInNeighbor[i], faceIndex, startIndex + int( i / 2 )} );
2252 
2253  int vertexIndex = face.at( ( splitVertexPos + i ) % faceSize );
2254  if ( mVertexToFace.at( vertexIndex ) == faceIndex )
2255  changes.mVerticesToFaceChanges.append( {vertexIndex, faceIndex, startIndex + int( i / 2 )} );
2256  }
2257 
2258  applyChanges( changes );
2259 
2260  return changes;
2261 }
2262 
2263 
2265 {
2266  Changes changes;
2267  changes.mVerticesToAdd.append( vertex );
2268  changes.mVertexToFaceToAdd.append( -1 );
2269 
2270  mMesh->vertices.append( vertex );
2271  mVertexToFace.append( -1 );
2272  changes.mAddedFacesFirstIndex = mMesh->faceCount();
2273 
2274  const QgsMeshFace includingFace = mMesh->face( includingFaceIndex );
2275  const FaceNeighbors includingFaceNeighborhood = mFacesNeighborhood.at( includingFaceIndex );
2276  int includingFaceSize = includingFace.count();
2277 
2278  for ( int i = 0; i < includingFaceSize; ++i )
2279  {
2280  // add a new face
2281  QgsMeshFace face( 3 );
2282  face[0] = mMesh->vertexCount() - 1;
2283  face[1] = includingFace.at( i );
2284  face[2] = includingFace.at( ( i + 1 ) % includingFaceSize );
2285  mMesh->faces.append( face );
2286  changes.mFacesToAdd.append( face );
2287 
2288  int currentVertexIndex = includingFace.at( i );
2289  if ( mVertexToFace.at( currentVertexIndex ) == includingFaceIndex )
2290  {
2291  int newFaceIndex = mMesh->faceCount() - 1;
2292  mVertexToFace[currentVertexIndex] = newFaceIndex;
2293  changes.mVerticesToFaceChanges.append( {currentVertexIndex, includingFaceIndex, newFaceIndex} );
2294  }
2295 
2296  int includingFaceNeighbor = includingFaceNeighborhood.at( i );
2297  FaceNeighbors neighbors(
2298  {
2299  changes.mAddedFacesFirstIndex + ( i + includingFaceSize - 1 ) % includingFaceSize,
2300  includingFaceNeighbor,
2301  changes.mAddedFacesFirstIndex + ( i + includingFaceSize + 1 ) % includingFaceSize
2302  } );
2303  mFacesNeighborhood.append( neighbors );
2304  changes.mFacesNeighborhoodToAdd.append( neighbors );
2305 
2306  if ( includingFaceNeighbor != -1 )
2307  {
2308  int indexInNeighbor = vertexPositionInFace( *mMesh, includingFace.at( ( i + 1 ) % includingFaceSize ), includingFaceNeighbor );
2309  int oldValue = mFacesNeighborhood[includingFaceNeighbor][indexInNeighbor];
2310  mFacesNeighborhood[includingFaceNeighbor][indexInNeighbor] = changes.mAddedFacesFirstIndex + i;
2311  changes.mNeighborhoodChanges.append( {includingFaceNeighbor, indexInNeighbor, oldValue, changes.mAddedFacesFirstIndex + i} );
2312  }
2313  }
2314 
2315  changes.mFacesToRemove.append( includingFace );
2316  changes.mFaceIndexesToRemove.append( includingFaceIndex );
2317  changes.mFacesNeighborhoodToRemove.append( includingFaceNeighborhood );
2318 
2319  mFacesNeighborhood[includingFaceIndex] = FaceNeighbors();
2320  mMesh->faces[includingFaceIndex] = QgsMeshFace();
2321  mVertexToFace[mVertexToFace.count() - 1] = mMesh->faceCount() - 1;
2322  changes.mVertexToFaceToAdd[changes.mVertexToFaceToAdd.count() - 1] = mMesh->faceCount() - 1 ;
2323 
2324  return changes;
2325 }
2326 
2328 {
2329  const QgsMeshFace face1 = mMesh->face( faceIndex );
2330 
2331  Changes changes;
2332  changes.mVerticesToAdd.append( vertexToInsert );
2333  changes.mAddedFacesFirstIndex = mMesh->faceCount();
2334 
2335  // triangulate the first face
2336  int newVertexPositionInFace1 = position + 1;
2337 
2338  auto triangulate = [this, &changes]( int removedFaceIndex, const QgsMeshVertex & newVertex, int newVertexPosition, QVector<int> &edgeFacesIndexes )->bool
2339  {
2340  const QgsMeshFace &initialFace = mMesh->face( removedFaceIndex );
2341  changes.mFacesToRemove.append( initialFace );
2342  changes.mFaceIndexesToRemove.append( removedFaceIndex );
2343  changes.mFacesNeighborhoodToRemove.append( mFacesNeighborhood.at( removedFaceIndex ) );
2344  const int addedVertexIndex = mMesh->vertexCount();
2345 
2346  int faceStartGlobalIndex = mMesh->faceCount() + changes.mFacesToAdd.count();
2347  int localStartIndex = changes.mFacesToAdd.count();
2348 
2349  QVector<int> newBoundary = initialFace;
2350  newBoundary.insert( newVertexPosition, addedVertexIndex );
2351 
2352  try
2353  {
2354  QHash<p2t::Point *, int> mapPoly2TriPointToVertex;
2355  std::vector<p2t::Point *> faceToFill( newBoundary.count() );
2356  for ( int i = 0; i < newBoundary.count(); ++i )
2357  {
2358  QgsMeshVertex vert;
2359 
2360  if ( newBoundary.at( i ) == addedVertexIndex )
2361  vert = newVertex;
2362  else
2363  vert = mMesh->vertex( newBoundary.at( i ) );
2364 
2365  faceToFill[i] = new p2t::Point( vert.x(), vert.y() );
2366  mapPoly2TriPointToVertex.insert( faceToFill[i], newBoundary.at( i ) );
2367  }
2368 
2369  std::unique_ptr<p2t::CDT> cdt( new p2t::CDT( faceToFill ) );
2370  cdt->Triangulate();
2371  std::vector<p2t::Triangle *> triangles = cdt->GetTriangles();
2372  QVector<QgsMeshFace> newFaces( triangles.size() );
2373  for ( size_t i = 0; i < triangles.size(); ++i )
2374  {
2375  QgsMeshFace &face = newFaces[i];
2376  face.resize( 3 );
2377  for ( int j = 0; j < 3; j++ )
2378  {
2379  int vertInd = mapPoly2TriPointToVertex.value( triangles.at( i )->GetPoint( j ), -1 );
2380  if ( vertInd == -1 )
2381  throw std::exception();
2382  face[j] = vertInd;
2383  }
2384  }
2385 
2386  QgsMeshEditingError error;
2387  QgsTopologicalMesh::TopologicalFaces topologicalFaces = createNewTopologicalFaces( newFaces, false, error );
2389  throw std::exception();
2390 
2391  changes.mFacesToAdd.append( topologicalFaces.meshFaces() );
2392  changes.mFacesNeighborhoodToAdd.append( topologicalFaces.mFacesNeighborhood );
2393 
2394  // vertices to face changes
2395  const QList<int> &verticesToFaceToChange = topologicalFaces.mVerticesToFace.keys();
2396  for ( const int vtc : verticesToFaceToChange )
2397  if ( vtc != addedVertexIndex && mVertexToFace.at( vtc ) == removedFaceIndex )
2398  changes.mVerticesToFaceChanges.append(
2399  {
2400  vtc,
2401  removedFaceIndex,
2402  topologicalFaces.vertexToFace( vtc ) + faceStartGlobalIndex
2403  } );
2404 
2405  // reindex neighborhood for new faces
2406  for ( int i = 0; i < topologicalFaces.mFaces.count(); ++i )
2407  {
2408  FaceNeighbors &faceNeighbors = changes.mFacesNeighborhoodToAdd[i + localStartIndex];
2409  for ( int n = 0; n < faceNeighbors.count(); ++n )
2410  {
2411  if ( faceNeighbors.at( n ) != -1 )
2412  faceNeighbors[n] += faceStartGlobalIndex; //reindex internal neighborhood
2413  }
2414  }
2415 
2416  edgeFacesIndexes.resize( 2 );
2417  // link neighborhood for boundaries of each side
2418  for ( int i = 0 ; i < newBoundary.count(); ++i )
2419  {
2420  int vertexIndex = newBoundary.at( i );
2421  QgsMeshVertexCirculator circulator = QgsMeshVertexCirculator( topologicalFaces, vertexIndex );
2422  circulator.goBoundaryClockwise();
2423  int newFaceBoundaryLocalIndex = localStartIndex + circulator.currentFaceIndex();
2424 
2425  int newFaceBoundaryIndexInMesh = faceStartGlobalIndex;
2426 
2427  int meshFaceBoundaryIndex;
2428  if ( i == newVertexPosition )
2429  {
2430  meshFaceBoundaryIndex = -1; //face that are on the opposite side of the edge, filled later
2431  edgeFacesIndexes[0] = newFaceBoundaryLocalIndex;
2432  }
2433  else if ( i == ( newVertexPosition + newBoundary.count() - 1 ) % newBoundary.count() )
2434  {
2435  meshFaceBoundaryIndex = -1; //face that are on the opposite side of the edge, filled later
2436  edgeFacesIndexes[1] = newFaceBoundaryLocalIndex;
2437  }
2438  else
2439  meshFaceBoundaryIndex = mFacesNeighborhood.at( removedFaceIndex ).at( vertexPositionInFace( vertexIndex, initialFace ) );
2440 
2441  const QgsMeshFace &newFace = circulator.currentFace();
2442  int positionInNewFaces = vertexPositionInFace( vertexIndex, newFace );
2443 
2444  if ( meshFaceBoundaryIndex != -1 )
2445  {
2446  const QgsMeshFace meshFace = mMesh->face( meshFaceBoundaryIndex );
2447  int positionInMeshFaceBoundary = vertexPositionInFace( *mMesh, vertexIndex, meshFaceBoundaryIndex );
2448  positionInMeshFaceBoundary = ( positionInMeshFaceBoundary - 1 + meshFace.count() ) % meshFace.count(); //take the position just before
2449 
2450  changes.mNeighborhoodChanges.append( {meshFaceBoundaryIndex,
2451  positionInMeshFaceBoundary,
2452  removedFaceIndex,
2453  newFaceBoundaryIndexInMesh +
2454  circulator.currentFaceIndex()} );
2455  }
2456 
2457  changes.mFacesNeighborhoodToAdd[newFaceBoundaryLocalIndex][positionInNewFaces] = meshFaceBoundaryIndex;
2458  }
2459 
2460  qDeleteAll( faceToFill );
2461  }
2462  catch ( ... )
2463  {
2464  return false;
2465  }
2466 
2467  return true;
2468  };
2469 
2470  QVector<int> edgeFacesIndexes;
2471  if ( !triangulate( faceIndex, vertexToInsert, newVertexPositionInFace1, edgeFacesIndexes ) )
2472  return Changes();
2473 
2474  changes.mVertexToFaceToAdd.append( edgeFacesIndexes.at( 0 ) + changes.mAddedFacesFirstIndex );
2475 
2476  int addedVertexIndex = mMesh->vertexCount();
2477 
2478  //if exist, triangulate the second face if exists
2479  int face2Index = mFacesNeighborhood.at( faceIndex ).at( position );
2480  if ( face2Index != -1 )
2481  {
2482  const QgsMeshFace &face2 = mMesh->face( face2Index );
2483  int vertexPositionInFace2 = vertexPositionInFace( face1.at( position ), face2 );
2484  QVector<int> edgeFacesIndexesFace2;
2485  if ( !triangulate( face2Index, vertexToInsert, vertexPositionInFace2, edgeFacesIndexesFace2 ) )
2486  return Changes();
2487 
2488  //link neighborhood with other side
2489  const QgsMeshFace &firstFaceSide1 = changes.mFacesToAdd.at( edgeFacesIndexes.at( 0 ) );
2490  int pos1InFaceSide1 = vertexPositionInFace( addedVertexIndex, firstFaceSide1 );
2491 
2492  const QgsMeshFace &secondFaceSide1 = changes.mFacesToAdd.at( edgeFacesIndexes.at( 1 ) );
2493  int pos2InFaceSide1 = vertexPositionInFace( addedVertexIndex, secondFaceSide1 );
2494  pos2InFaceSide1 = ( pos2InFaceSide1 + secondFaceSide1.size() - 1 ) % secondFaceSide1.size();
2495 
2496  const QgsMeshFace &firstFaceSide2 = changes.mFacesToAdd.at( edgeFacesIndexesFace2.at( 0 ) );
2497  int pos1InFaceSide2 = vertexPositionInFace( addedVertexIndex, firstFaceSide2 );
2498 
2499  const QgsMeshFace &secondFaceSide2 = changes.mFacesToAdd.at( edgeFacesIndexesFace2.at( 1 ) );
2500  int pos2InFaceSide2 = vertexPositionInFace( addedVertexIndex, secondFaceSide2 );
2501  pos2InFaceSide2 = ( pos2InFaceSide2 + secondFaceSide1.size() - 1 ) % secondFaceSide1.size();
2502 
2503  changes.mFacesNeighborhoodToAdd[edgeFacesIndexes.at( 0 )][pos1InFaceSide1] = edgeFacesIndexesFace2.at( 1 ) + changes.mAddedFacesFirstIndex;
2504  changes.mFacesNeighborhoodToAdd[edgeFacesIndexes.at( 1 )][pos2InFaceSide1] = edgeFacesIndexesFace2.at( 0 ) + changes.mAddedFacesFirstIndex;
2505  changes.mFacesNeighborhoodToAdd[edgeFacesIndexesFace2.at( 0 )][pos1InFaceSide2] = edgeFacesIndexes.at( 1 ) + changes.mAddedFacesFirstIndex;
2506  changes.mFacesNeighborhoodToAdd[edgeFacesIndexesFace2.at( 1 )][pos2InFaceSide2] = edgeFacesIndexes.at( 0 ) + changes.mAddedFacesFirstIndex;
2507  }
2508 
2509  applyChanges( changes );
2510  return changes;
2511 }
2512 
2513 QgsTopologicalMesh::Changes QgsTopologicalMesh::changeZValue( const QList<int> &verticesIndexes, const QList<double> &newValues )
2514 {
2515  Q_ASSERT( verticesIndexes.count() == newValues.count() );
2516  Changes changes;
2517  changes.mChangeCoordinateVerticesIndexes.reserve( verticesIndexes.count() );
2518  changes.mNewZValues.reserve( verticesIndexes.count() );
2519  changes.mOldZValues.reserve( verticesIndexes.count() );
2520  for ( int i = 0; i < verticesIndexes.count(); ++i )
2521  {
2522  changes.mChangeCoordinateVerticesIndexes.append( verticesIndexes.at( i ) );
2523  changes.mOldZValues.append( mMesh->vertices.at( verticesIndexes.at( i ) ).z() );
2524  changes.mNewZValues.append( newValues.at( i ) );
2525  }
2526 
2527  applyChanges( changes );
2528 
2529  return changes;
2530 }
2531 
2532 QgsTopologicalMesh::Changes QgsTopologicalMesh::changeXYValue( const QList<int> &verticesIndexes, const QList<QgsPointXY> &newValues )
2533 {
2534  Q_ASSERT( verticesIndexes.count() == newValues.count() );
2535  Changes changes;
2536  changes.mChangeCoordinateVerticesIndexes.reserve( verticesIndexes.count() );
2537  changes.mNewXYValues.reserve( verticesIndexes.count() );
2538  changes.mOldXYValues.reserve( verticesIndexes.count() );
2539  QSet<int> concernedFace;
2540  for ( int i = 0; i < verticesIndexes.count(); ++i )
2541  {
2542  changes.mChangeCoordinateVerticesIndexes.append( verticesIndexes.at( i ) );
2543  changes.mOldXYValues.append( mMesh->vertices.at( verticesIndexes.at( i ) ) );
2544  changes.mNewXYValues.append( newValues.at( i ) );
2545  concernedFace.unite( qgis::listToSet( facesAroundVertex( verticesIndexes.at( i ) ) ) );
2546  }
2547 
2548  changes.mNativeFacesIndexesGeometryChanged = concernedFace.values();
2549 
2550  applyChanges( changes );
2551 
2552  return changes;
2553 }
@ InvalidFace
An error occurs due to an invalid face (for example, vertex indexes are unordered)
@ UniqueSharedVertex
A least two faces share only one vertices.
@ ManifoldFace
ManifoldFace.
@ InvalidVertex
An error occurs due to an invalid vertex (for example, vertex index is out of range the available ver...
@ FlatFace
A flat face is present.
static bool segmentIntersection(const QgsPoint &p1, const QgsPoint &p2, const QgsPoint &q1, const QgsPoint &q2, QgsPoint &intersectionPoint, bool &isIntersection, double tolerance=1e-8, bool acceptImproperIntersection=false) SIP_HOLDGIL
Compute the intersection between two segments.
Class that represents an error during mesh editing.
Definition: qgsmesheditor.h:43
Qgis::MeshEditingErrorType errorType
Definition: qgsmesheditor.h:52
Convenient class that turn around a vertex and provide information about faces and vertices.
bool isValid() const
Returns whether the vertex circulator is valid.
int turnClockwise() const
Turns counter clockwise around the vertex and returns the new current face, -1 if the circulator pass...
bool goBoundaryCounterClockwise() const
Sets the circulator on the boundary face turning counter clockwise, return false is there isn't bound...
int oppositeVertexCounterClockwise() const
Returns the opposite vertex of the current face and on the edge on the side turning counter clockwise...
int turnCounterClockwise() const
Turns counter clockwise around the vertex and returns the new current face, -1 if the circulator pass...
int currentFaceIndex() const
Returns the current face index, -1 if the circulator has passed a boundary or circulator is invalid.
bool goBoundaryClockwise() const
Sets the circulator on the boundary face turning clockwise, return false is there isn't boundary face...
QgsMeshFace currentFace() const
Returns the current face, empty face if the circulator pass a boundary or circulator is invalid.
QgsMeshVertexCirculator(const QgsTopologicalMesh &topologicalMesh, int vertexIndex)
Constructor with topologicalMesh and vertexIndex.
int oppositeVertexClockwise() const
Returns the opposite vertex of the current face and on the edge on the side turning clockwise.
int degree() const
Returns the degree of the vertex, that is the count of other vertices linked.
QList< int > facesAround() const
Returns all the faces indexes around the vertex.
static void logMessage(const QString &message, const QString &tag=QString(), Qgis::MessageLevel level=Qgis::MessageLevel::Warning, bool notifyUser=true)
Adds a message to the log instance (and creates it if necessary).
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
bool isEmpty() const override SIP_HOLDGIL
Returns true if the geometry is empty.
Definition: qgspoint.cpp:767
Q_GADGET double x
Definition: qgspoint.h:52
double y
Definition: qgspoint.h:53
Class that contains topological differences between two states of a topological mesh,...
QList< int > mChangeCoordinateVerticesIndexes
void clearChanges()
Clears all changes.
QVector< FaceNeighbors > mFacesNeighborhoodToRemove
QVector< QgsMeshFace > removedFaces() const
Returns the faces that are removed with this changes.
QList< QgsPointXY > mNewXYValues
QList< QgsMeshVertex > mRemovedVertices
QVector< QgsMeshVertex > addedVertices() const
Returns the added vertices with this changes.
QList< std::array< int, 4 > > mNeighborhoodChanges
bool isEmpty() const
Returns whether changes are empty, that there is nothing to change.
QList< int > changedCoordinatesVerticesIndexes() const
Returns the indexes of vertices that have changed coordinates.
QList< int > mNativeFacesIndexesGeometryChanged
QVector< QgsMeshFace > mFacesToAdd
QList< int > removedFaceIndexes() const
Returns the indexes of the faces that are removed with this changes.
QVector< FaceNeighbors > mFacesNeighborhoodToAdd
QList< std::array< int, 3 > > mVerticesToFaceChanges
QList< QgsPointXY > mOldXYValues
QList< double > newVerticesZValues() const
Returns the new Z values of vertices that have changed their coordinates.
QVector< QgsMeshVertex > mVerticesToAdd
QVector< QgsMeshFace > addedFaces() const
Returns the face that are added with this changes.
QList< QgsPointXY > oldVerticesXYValues() const
Returns the old (X,Y) values of vertices that have changed their coordinates.
QList< QgsPointXY > newVerticesXYValues() const
Returns the new (X,Y) values of vertices that have changed their coordinates.
QVector< QgsMeshFace > mFacesToRemove
QList< int > nativeFacesIndexesGeometryChanged() const
Returns a list of the native face indexes that have a geometry changed.
QList< int > verticesToRemoveIndexes() const
Returns the indexes of vertices to remove.
Class that contains independent faces an topological information about this faces.
int vertexToFace(int vertexIndex) const
Returns a face linked to the vertices with index vertexIndex.
QVector< FaceNeighbors > facesNeighborhood() const
Returns the face neighborhood of the faces, indexing is local.
void clear()
Clears all data contained in the instance.
QVector< QgsMeshFace > meshFaces() const
Returns faces.
Class that wraps a QgsMesh to ensure the consistency of the mesh during editing and help to access to...
Changes changeZValue(const QList< int > &verticesIndexes, const QList< double > &newValues)
Changes the Z values of the vertices with indexes in vertices indexes with the values in newValues.
static QgsTopologicalMesh createTopologicalMesh(QgsMesh *mesh, int maxVerticesPerFace, QgsMeshEditingError &error)
Creates a topologicaly consistent mesh with mesh, this static method modifies mesh to be topological ...
bool isVertexFree(int vertexIndex) const
Returns whether the vertex is a free vertex.
static QgsMeshEditingError counterClockwiseFaces(QgsMeshFace &face, QgsMesh *mesh)
Checks the topology of the face and sets it counter clockwise if necessary.
Changes removeVertexFillHole(int vertexIndex)
Removes the vertex with index vertexIndex.
static QgsMeshEditingError checkTopology(const QgsMesh &mesh, int maxVerticesPerFace)
Checks the topology of the mesh mesh, if error occurs, this mesh can't be edited.
friend class QgsMeshVertexCirculator
void applyChanges(const Changes &changes)
Applies the changes.
int firstFaceLinked(int vertexIndex) const
Returns the index of the first face linked, returns -1 if it is a free vertex or out of range index.
QgsMeshEditingError facesCanBeRemoved(const QList< int > facesIndexes)
Returns whether faces with index in faceIndexes can be removed/ The method an error object with type ...
QgsMeshEditingError checkConsistency() const
Checks the consistency of the topological mesh and return false if there is a consistency issue.
Changes removeVertices(const QList< int > &vertices)
Removes all the vertices with index in the list vertices If vertices in linked with faces,...
Changes changeXYValue(const QList< int > &verticesIndexes, const QList< QgsPointXY > &newValues)
Changes the (X,Y) values of the vertices with indexes in vertices indexes with the values in newValue...
void reindex()
Reindexes faces and vertices, after this operation, the topological mesh can't be edited anymore and ...
QVector< int > neighborsOfFace(int faceIndex) const
Returns the indexes of neighbor faces of the face with index faceIndex.
QgsMeshEditingError facesCanBeAdded(const TopologicalFaces &topologicalFaces) const
Returns whether the faces can be added to the mesh.
bool renumber()
Renumbers the indexes of vertices and faces using the Reverse CutHill McKee Algorithm.
Changes flipEdge(int vertexIndex1, int vertexIndex2)
Flips edge (vertexIndex1, vertexIndex2) The method returns a instance of the class QgsTopologicalMesh...
QVector< int > FaceNeighbors
void reverseChanges(const Changes &changes)
Reverses the changes.
Changes addFaces(const TopologicalFaces &topologicFaces)
Adds faces topologicFaces to the topologic mesh.
Changes merge(int vertexIndex1, int vertexIndex2)
Merges faces separated by vertices with indexes vertexIndex1 and vertexIndex2 The method returns a in...
QList< int > freeVerticesIndexes() const
Returns a list of vertices are not linked to any faces.
bool edgeCanBeFlipped(int vertexIndex1, int vertexIndex2) const
Returns true if the edge can be flipped (only available for edge shared by two faces with 3 vertices)
Changes addVertexInFace(int faceIndex, const QgsMeshVertex &vertex)
Adds a vertex in the face with index faceIndex.
Changes removeFaces(const QList< int > facesIndexes)
Removes faces with index in faceIndexes.
bool canBeMerged(int vertexIndex1, int vertexIndex2) const
Returns true if faces separated by vertices with indexes vertexIndex1 and vertexIndex2 can be merged.
QList< int > facesAroundVertex(int vertexIndex) const
Returns the indexes of faces that are around the vertex with index vertexIndex.
bool canBeSplit(int faceIndex) const
Returns true if face with index faceIndex can be split.
Changes addFreeVertex(const QgsMeshVertex &vertex)
Adds a free vertex in the face, that is a vertex tha tis not included or linked with any faces.
Changes insertVertexInFacesEdge(int faceIndex, int position, const QgsMeshVertex &vertex)
Inserts a vertex in the edge of face with index faceIndex at position .
QgsMesh * mesh() const
Returns a pointer to the wrapped mesh.
bool isVertexOnBoundary(int vertexIndex) const
Returns whether the vertex is on a boundary.
Changes splitFace(int faceIndex)
Splits face with index faceIndex The method returns a instance of the class QgsTopologicalMesh::Chang...
static TopologicalFaces createNewTopologicalFaces(const QVector< QgsMeshFace > &faces, bool uniqueSharedVertexAllowed, QgsMeshEditingError &error)
Creates new topological faces that are not yet included in the mesh.
QgsMeshVertexCirculator vertexCirculator(int vertexIndex) const
Returns a vertex circulator linked to this mesh around the vertex with index vertexIndex.
A class to represent a vector.
Definition: qgsvector.h:30
double angle() const SIP_HOLDGIL
Returns the angle of the vector in radians.
Definition: qgsvector.h:172
double ANALYSIS_EXPORT angle(QgsPoint *p1, QgsPoint *p2, QgsPoint *p3, QgsPoint *p4)
Calculates the angle between two segments (in 2 dimension, z-values are ignored)
Definition: MathUtils.cpp:786
QVector< int > QgsMeshFace
List of vertex indexes.
QgsPoint QgsMeshVertex
xyz coords of vertex
Mesh - vertices, edges and faces.
int vertexCount() const
Returns number of vertices.
QVector< QgsMeshVertex > vertices
QgsMeshFace face(int index) const
Returns a face at the index.
QVector< QgsMeshFace > faces
int faceCount() const
Returns number of faces.
QgsMeshVertex vertex(int index) const
Returns a vertex at the index.