QGIS API Documentation 3.34.0-Prizren (ffbdd678812)
Loading...
Searching...
No Matches
qgsinternalgeometryengine.cpp
Go to the documentation of this file.
1/***************************************************************************
2 qgsinternalgeometryengine.cpp - QgsInternalGeometryEngine
3
4 ---------------------
5 begin : 13.1.2016
6 copyright : (C) 2016 by Matthias Kuhn
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
18
19#include "qgslinestring.h"
20#include "qgsmultipolygon.h"
21#include "qgspolygon.h"
22#include "qgsmulticurve.h"
23#include "qgscircularstring.h"
24#include "qgsgeometry.h"
25#include "qgsgeometryutils.h"
26#include "qgslinesegment.h"
27#include "qgscircle.h"
28#include "qgstessellator.h"
29#include "qgsfeedback.h"
30#include "qgsgeometryengine.h"
31#include "qgsmultilinestring.h"
32#include <QTransform>
33#include <functional>
34#include <memory>
35#include <queue>
36#include <random>
37
39 : mGeometry( geometry.constGet() )
40{
41
42}
43
45{
46 return mLastError;
47}
48
49
50enum class Direction
51{
52 Up,
53 Right,
54 Down,
55 Left,
56 None
57};
58
64Direction getEdgeDirection( const QgsPoint &p1, const QgsPoint &p2, double maxDev )
65{
66 double dx = p2.x() - p1.x();
67 double dy = p2.y() - p1.y();
68 if ( ( dx == 0.0 ) && ( dy == 0.0 ) )
69 return Direction::None;
70 if ( fabs( dx ) >= fabs( dy ) )
71 {
72 double dev = fabs( dy ) / fabs( dx );
73 if ( dev > maxDev )
74 return Direction::None;
75 return dx > 0.0 ? Direction::Right : Direction::Left;
76 }
77 else
78 {
79 double dev = fabs( dx ) / fabs( dy );
80 if ( dev > maxDev )
81 return Direction::None;
82 return dy > 0.0 ? Direction::Up : Direction::Down;
83 }
84}
85
91std::pair<bool, std::array<Direction, 4>> getEdgeDirections( const QgsPolygon *g, double maxDev )
92{
93 std::pair<bool, std::array<Direction, 4>> ret = { false, { Direction::None, Direction::None, Direction::None, Direction::None } };
94 // The polygon might start in the middle of a side. Hence, we need a fifth
95 // direction to record the beginning of the side when we went around the
96 // polygon.
97 std::array<Direction, 5> dirs;
98
99 int idx = 0;
101 QgsAbstractGeometry::vertex_iterator current = previous;
102 ++current;
104 while ( current != end )
105 {
106 Direction dir = getEdgeDirection( *previous, *current, maxDev );
107 if ( dir == Direction::None )
108 return ret;
109 if ( idx == 0 )
110 {
111 dirs[0] = dir;
112 ++idx;
113 }
114 else if ( dir != dirs[idx - 1] )
115 {
116 if ( idx == 5 )
117 return ret;
118 dirs[idx] = dir;
119 ++idx;
120 }
121 previous = current;
122 ++current;
123 }
124 ret.first = ( idx == 5 ) ? ( dirs[0] == dirs[4] ) : ( idx == 4 );
125 std::copy( dirs.begin(), dirs.begin() + 4, ret.second.begin() );
126 return ret;
127}
128
129bool matchesOrientation( std::array<Direction, 4> dirs, std::array<Direction, 4> oriented )
130{
131 int idx = std::find( oriented.begin(), oriented.end(), dirs[0] ) - oriented.begin();
132 for ( int i = 1; i < 4; ++i )
133 {
134 if ( dirs[i] != oriented[( idx + i ) % 4] )
135 return false;
136 }
137 return true;
138}
139
143bool isClockwise( std::array<Direction, 4> dirs )
144{
145 const std::array<Direction, 4> cwdirs = { Direction::Up, Direction::Right, Direction::Down, Direction::Left };
146 return matchesOrientation( dirs, cwdirs );
147}
148
153bool isCounterClockwise( std::array<Direction, 4> dirs )
154{
155 const std::array<Direction, 4> ccwdirs = { Direction::Right, Direction::Up, Direction::Left, Direction::Down };
156 return matchesOrientation( dirs, ccwdirs );
157}
158
159
160bool QgsInternalGeometryEngine::isAxisParallelRectangle( double maximumDeviation, bool simpleRectanglesOnly ) const
161{
163 return false;
164
165 const QgsPolygon *polygon = qgsgeometry_cast< const QgsPolygon * >( mGeometry );
166 if ( !polygon->exteriorRing() || polygon->numInteriorRings() > 0 )
167 return false;
168
169 const int vertexCount = polygon->exteriorRing()->numPoints();
170 if ( vertexCount < 4 )
171 return false;
172 else if ( simpleRectanglesOnly && ( vertexCount != 5 || !polygon->exteriorRing()->isClosed() ) )
173 return false;
174
175 bool found4Dirs;
176 std::array<Direction, 4> dirs;
177 std::tie( found4Dirs, dirs ) = getEdgeDirections( polygon, std::tan( maximumDeviation * M_PI / 180 ) );
178
179 return found4Dirs && ( isCounterClockwise( dirs ) || isClockwise( dirs ) );
180}
181
182/***************************************************************************
183 * This class is considered CRITICAL and any change MUST be accompanied with
184 * full unit tests.
185 * See details in QEP #17
186 ****************************************************************************/
187
189{
190 mLastError.clear();
191 QVector<QgsLineString *> linesToProcess;
192
193 const QgsMultiCurve *multiCurve = qgsgeometry_cast< const QgsMultiCurve * >( mGeometry );
194 if ( multiCurve )
195 {
196 linesToProcess.reserve( multiCurve->partCount() );
197 for ( int i = 0; i < multiCurve->partCount(); ++i )
198 {
199 linesToProcess << static_cast<QgsLineString *>( multiCurve->geometryN( i )->clone() );
200 }
201 }
202
203 const QgsCurve *curve = qgsgeometry_cast< const QgsCurve * >( mGeometry );
204 if ( curve )
205 {
206 linesToProcess << static_cast<QgsLineString *>( curve->segmentize() );
207 }
208
209 std::unique_ptr<QgsMultiPolygon> multipolygon( linesToProcess.size() > 1 ? new QgsMultiPolygon() : nullptr );
210 QgsPolygon *polygon = nullptr;
211
212 if ( !linesToProcess.empty() )
213 {
214 std::unique_ptr< QgsLineString > secondline;
215 for ( QgsLineString *line : std::as_const( linesToProcess ) )
216 {
217 QTransform transform = QTransform::fromTranslate( x, y );
218
219 secondline.reset( line->reversed() );
220 secondline->transform( transform );
221
222 line->append( secondline.get() );
223 line->addVertex( line->pointN( 0 ) );
224
225 polygon = new QgsPolygon();
226 polygon->setExteriorRing( line );
227
228 if ( multipolygon )
229 multipolygon->addGeometry( polygon );
230 }
231
232 if ( multipolygon )
233 return QgsGeometry( multipolygon.release() );
234 else
235 return QgsGeometry( polygon );
236 }
237
238 return QgsGeometry();
239}
240
241
242
243// polylabel implementation
244// ported from the original JavaScript implementation developed by Vladimir Agafonkin
245// originally licensed under the ISC License
246
248class Cell
249{
250 public:
251 Cell( double x, double y, double h, const QgsPolygon *polygon )
252 : x( x )
253 , y( y )
254 , h( h )
255 , d( polygon->pointDistanceToBoundary( x, y ) )
256 , max( d + h * M_SQRT2 )
257 {}
258
260 double x;
262 double y;
264 double h;
266 double d;
268 double max;
269};
270
271struct GreaterThanByMax
272{
273 bool operator()( const Cell *lhs, const Cell *rhs ) const
274 {
275 return rhs->max > lhs->max;
276 }
277};
278
279Cell *getCentroidCell( const QgsPolygon *polygon )
280{
281 double area = 0;
282 double x = 0;
283 double y = 0;
284
285 const QgsLineString *exterior = static_cast< const QgsLineString *>( polygon->exteriorRing() );
286 int len = exterior->numPoints() - 1; //assume closed
287 for ( int i = 0, j = len - 1; i < len; j = i++ )
288 {
289 double aX = exterior->xAt( i );
290 double aY = exterior->yAt( i );
291 double bX = exterior->xAt( j );
292 double bY = exterior->yAt( j );
293 double f = aX * bY - bX * aY;
294 x += ( aX + bX ) * f;
295 y += ( aY + bY ) * f;
296 area += f * 3;
297 }
298 if ( area == 0 )
299 return new Cell( exterior->xAt( 0 ), exterior->yAt( 0 ), 0, polygon );
300 else
301 return new Cell( x / area, y / area, 0.0, polygon );
302}
303
304QgsPoint surfacePoleOfInaccessibility( const QgsSurface *surface, double precision, double &distanceFromBoundary )
305{
306 std::unique_ptr< QgsPolygon > segmentizedPoly;
307 const QgsPolygon *polygon = qgsgeometry_cast< const QgsPolygon * >( surface );
308 if ( !polygon )
309 {
310 segmentizedPoly.reset( static_cast< QgsPolygon *>( surface->segmentize() ) );
311 polygon = segmentizedPoly.get();
312 }
313
314 // start with the bounding box
315 QgsRectangle bounds = polygon->boundingBox();
316
317 // initial parameters
318 double cellSize = std::min( bounds.width(), bounds.height() );
319
320 if ( qgsDoubleNear( cellSize, 0.0 ) )
321 return QgsPoint( bounds.xMinimum(), bounds.yMinimum() );
322
323 double h = cellSize / 2.0;
324 std::priority_queue< Cell *, std::vector<Cell *>, GreaterThanByMax > cellQueue;
325
326 // cover polygon with initial cells
327 for ( double x = bounds.xMinimum(); x < bounds.xMaximum(); x += cellSize )
328 {
329 for ( double y = bounds.yMinimum(); y < bounds.yMaximum(); y += cellSize )
330 {
331 cellQueue.push( new Cell( x + h, y + h, h, polygon ) );
332 }
333 }
334
335 // take centroid as the first best guess
336 std::unique_ptr< Cell > bestCell( getCentroidCell( polygon ) );
337
338 // special case for rectangular polygons
339 std::unique_ptr< Cell > bboxCell( new Cell( bounds.xMinimum() + bounds.width() / 2.0,
340 bounds.yMinimum() + bounds.height() / 2.0,
341 0, polygon ) );
342 if ( bboxCell->d > bestCell->d )
343 {
344 bestCell = std::move( bboxCell );
345 }
346
347 while ( !cellQueue.empty() )
348 {
349 // pick the most promising cell from the queue
350 std::unique_ptr< Cell > cell( cellQueue.top() );
351 cellQueue.pop();
352 Cell *currentCell = cell.get();
353
354 // update the best cell if we found a better one
355 if ( currentCell->d > bestCell->d )
356 {
357 bestCell = std::move( cell );
358 }
359
360 // do not drill down further if there's no chance of a better solution
361 if ( currentCell->max - bestCell->d <= precision )
362 continue;
363
364 // split the cell into four cells
365 h = currentCell->h / 2.0;
366 cellQueue.push( new Cell( currentCell->x - h, currentCell->y - h, h, polygon ) );
367 cellQueue.push( new Cell( currentCell->x + h, currentCell->y - h, h, polygon ) );
368 cellQueue.push( new Cell( currentCell->x - h, currentCell->y + h, h, polygon ) );
369 cellQueue.push( new Cell( currentCell->x + h, currentCell->y + h, h, polygon ) );
370 }
371
372 distanceFromBoundary = bestCell->d;
373
374 return QgsPoint( bestCell->x, bestCell->y );
375}
376
378
379QgsGeometry QgsInternalGeometryEngine::poleOfInaccessibility( double precision, double *distanceFromBoundary ) const
380{
381 mLastError.clear();
382 if ( distanceFromBoundary )
383 *distanceFromBoundary = std::numeric_limits<double>::max();
384
385 if ( !mGeometry || mGeometry->isEmpty() )
386 return QgsGeometry();
387
388 if ( precision <= 0 )
389 return QgsGeometry();
390
391 if ( const QgsGeometryCollection *gc = qgsgeometry_cast< const QgsGeometryCollection *>( mGeometry ) )
392 {
393 int numGeom = gc->numGeometries();
394 double maxDist = 0;
395 QgsPoint bestPoint;
396 for ( int i = 0; i < numGeom; ++i )
397 {
398 const QgsSurface *surface = qgsgeometry_cast< const QgsSurface * >( gc->geometryN( i ) );
399 if ( !surface )
400 continue;
401
402 double dist = std::numeric_limits<double>::max();
403 QgsPoint p = surfacePoleOfInaccessibility( surface, precision, dist );
404 if ( dist > maxDist )
405 {
406 maxDist = dist;
407 bestPoint = p;
408 }
409 }
410
411 if ( bestPoint.isEmpty() )
412 return QgsGeometry();
413
414 if ( distanceFromBoundary )
415 *distanceFromBoundary = maxDist;
416 return QgsGeometry( new QgsPoint( bestPoint ) );
417 }
418 else
419 {
420 const QgsSurface *surface = qgsgeometry_cast< const QgsSurface * >( mGeometry );
421 if ( !surface )
422 return QgsGeometry();
423
424 double dist = std::numeric_limits<double>::max();
425 QgsPoint p = surfacePoleOfInaccessibility( surface, precision, dist );
426 if ( p.isEmpty() )
427 return QgsGeometry();
428
429 if ( distanceFromBoundary )
430 *distanceFromBoundary = dist;
431 return QgsGeometry( new QgsPoint( p ) );
432 }
433}
434
435
436// helpers for orthogonalize
437// adapted from original code in potlatch/id osm editor
438
439bool dotProductWithinAngleTolerance( double dotProduct, double lowerThreshold, double upperThreshold )
440{
441 return lowerThreshold > std::fabs( dotProduct ) || std::fabs( dotProduct ) > upperThreshold;
442}
443
444double normalizedDotProduct( const QgsPoint &a, const QgsPoint &b, const QgsPoint &c )
445{
446 QgsVector p = a - b;
447 QgsVector q = c - b;
448
449 if ( p.length() > 0 )
450 p = p.normalized();
451 if ( q.length() > 0 )
452 q = q.normalized();
453
454 return p * q;
455}
456
457double squareness( QgsLineString *ring, double lowerThreshold, double upperThreshold )
458{
459 double sum = 0.0;
460
461 bool isClosed = ring->isClosed();
462 int numPoints = ring->numPoints();
463 QgsPoint a;
464 QgsPoint b;
465 QgsPoint c;
466
467 for ( int i = 0; i < numPoints; ++i )
468 {
469 if ( !isClosed && i == numPoints - 1 )
470 break;
471 else if ( !isClosed && i == 0 )
472 {
473 b = ring->pointN( 0 );
474 c = ring->pointN( 1 );
475 }
476 else
477 {
478 if ( i == 0 )
479 {
480 a = ring->pointN( numPoints - 1 );
481 b = ring->pointN( 0 );
482 }
483 if ( i == numPoints - 1 )
484 c = ring->pointN( 0 );
485 else
486 c = ring->pointN( i + 1 );
487
488 double dotProduct = normalizedDotProduct( a, b, c );
489 if ( !dotProductWithinAngleTolerance( dotProduct, lowerThreshold, upperThreshold ) )
490 continue;
491
492 sum += 2.0 * std::min( std::fabs( dotProduct - 1.0 ), std::min( std::fabs( dotProduct ), std::fabs( dotProduct + 1 ) ) );
493 }
494 a = b;
495 b = c;
496 }
497
498 return sum;
499}
500
501QgsVector calcMotion( const QgsPoint &a, const QgsPoint &b, const QgsPoint &c,
502 double lowerThreshold, double upperThreshold )
503{
504 QgsVector p = a - b;
505 QgsVector q = c - b;
506
507 if ( qgsDoubleNear( p.length(), 0.0 ) || qgsDoubleNear( q.length(), 0.0 ) )
508 return QgsVector( 0, 0 );
509
510 // 2.0 is a magic number from the original JOSM source code
511 double scale = 2.0 * std::min( p.length(), q.length() );
512
513 p = p.normalized();
514 q = q.normalized();
515 double dotProduct = p * q;
516
517 if ( !dotProductWithinAngleTolerance( dotProduct, lowerThreshold, upperThreshold ) )
518 {
519 return QgsVector( 0, 0 );
520 }
521
522 // wonderful nasty hack which has survived through JOSM -> id -> QGIS
523 // to deal with almost-straight segments (angle is closer to 180 than to 90/270).
524 if ( dotProduct < -M_SQRT1_2 )
525 dotProduct += 1.0;
526
527 QgsVector new_v = p + q;
528 if ( qgsDoubleNear( new_v.length(), 0.0 ) )
529 {
530 return QgsVector( 0, 0 );
531 }
532 // 0.1 magic number from JOSM implementation - think this is to limit each iterative step
533 return new_v.normalized() * ( 0.1 * dotProduct * scale );
534}
535
536QgsLineString *doOrthogonalize( QgsLineString *ring, int iterations, double tolerance, double lowerThreshold, double upperThreshold )
537{
538 double minScore = std::numeric_limits<double>::max();
539
540 bool isClosed = ring->isClosed();
541 int numPoints = ring->numPoints();
542
543 std::unique_ptr< QgsLineString > best( ring->clone() );
544
545 QVector< QgsVector > /* yep */ motions;
546 motions.reserve( numPoints );
547
548 for ( int it = 0; it < iterations; ++it )
549 {
550 motions.resize( 0 ); // avoid re-allocations
551
552 // first loop through an calculate all motions
553 QgsPoint a;
554 QgsPoint b;
555 QgsPoint c;
556 for ( int i = 0; i < numPoints; ++i )
557 {
558 if ( isClosed && i == numPoints - 1 )
559 motions << motions.at( 0 ); //closed ring, so last point follows first point motion
560 else if ( !isClosed && ( i == 0 || i == numPoints - 1 ) )
561 {
562 b = ring->pointN( 0 );
563 c = ring->pointN( 1 );
564 motions << QgsVector( 0, 0 ); //non closed line, leave first/last vertex untouched
565 }
566 else
567 {
568 if ( i == 0 )
569 {
570 a = ring->pointN( numPoints - 1 );
571 b = ring->pointN( 0 );
572 }
573 if ( i == numPoints - 1 )
574 c = ring->pointN( 0 );
575 else
576 c = ring->pointN( i + 1 );
577
578 motions << calcMotion( a, b, c, lowerThreshold, upperThreshold );
579 }
580 a = b;
581 b = c;
582 }
583
584 // then apply them
585 for ( int i = 0; i < numPoints; ++i )
586 {
587 ring->setXAt( i, ring->xAt( i ) + motions.at( i ).x() );
588 ring->setYAt( i, ring->yAt( i ) + motions.at( i ).y() );
589 }
590
591 double newScore = squareness( ring, lowerThreshold, upperThreshold );
592 if ( newScore < minScore )
593 {
594 best.reset( ring->clone() );
595 minScore = newScore;
596 }
597
598 if ( minScore < tolerance )
599 break;
600 }
601
602 delete ring;
603
604 return best.release();
605}
606
607
608QgsAbstractGeometry *orthogonalizeGeom( const QgsAbstractGeometry *geom, int maxIterations, double tolerance, double lowerThreshold, double upperThreshold )
609{
610 std::unique_ptr< QgsAbstractGeometry > segmentizedCopy;
611 if ( QgsWkbTypes::isCurvedType( geom->wkbType() ) )
612 {
613 segmentizedCopy.reset( geom->segmentize() );
614 geom = segmentizedCopy.get();
615 }
616
618 {
619 return doOrthogonalize( static_cast< QgsLineString * >( geom->clone() ),
620 maxIterations, tolerance, lowerThreshold, upperThreshold );
621 }
622 else
623 {
624 // polygon
625 const QgsPolygon *polygon = static_cast< const QgsPolygon * >( geom );
626 QgsPolygon *result = new QgsPolygon();
627
628 result->setExteriorRing( doOrthogonalize( static_cast< QgsLineString * >( polygon->exteriorRing()->clone() ),
629 maxIterations, tolerance, lowerThreshold, upperThreshold ) );
630 for ( int i = 0; i < polygon->numInteriorRings(); ++i )
631 {
632 result->addInteriorRing( doOrthogonalize( static_cast< QgsLineString * >( polygon->interiorRing( i )->clone() ),
633 maxIterations, tolerance, lowerThreshold, upperThreshold ) );
634 }
635
636 return result;
637 }
638}
639
640QgsGeometry QgsInternalGeometryEngine::orthogonalize( double tolerance, int maxIterations, double angleThreshold ) const
641{
642 mLastError.clear();
643 if ( !mGeometry || ( QgsWkbTypes::geometryType( mGeometry->wkbType() ) != Qgis::GeometryType::Line
645 {
646 return QgsGeometry();
647 }
648
649 double lowerThreshold = std::cos( ( 90 - angleThreshold ) * M_PI / 180.00 );
650 double upperThreshold = std::cos( angleThreshold * M_PI / 180.0 );
651
652 if ( const QgsGeometryCollection *gc = qgsgeometry_cast< const QgsGeometryCollection *>( mGeometry ) )
653 {
654 int numGeom = gc->numGeometries();
655 QVector< QgsAbstractGeometry * > geometryList;
656 geometryList.reserve( numGeom );
657 for ( int i = 0; i < numGeom; ++i )
658 {
659 geometryList << orthogonalizeGeom( gc->geometryN( i ), maxIterations, tolerance, lowerThreshold, upperThreshold );
660 }
661
662 QgsGeometry first = QgsGeometry( geometryList.takeAt( 0 ) );
663 for ( QgsAbstractGeometry *g : std::as_const( geometryList ) )
664 {
665 first.addPart( g );
666 }
667 return first;
668 }
669 else
670 {
671 return QgsGeometry( orthogonalizeGeom( mGeometry, maxIterations, tolerance, lowerThreshold, upperThreshold ) );
672 }
673}
674
675// if extraNodesPerSegment < 0, then use distance based mode
676QgsLineString *doDensify( const QgsLineString *ring, int extraNodesPerSegment = -1, double distance = 1 )
677{
678 QVector< double > outX;
679 QVector< double > outY;
680 QVector< double > outZ;
681 QVector< double > outM;
682 double multiplier = 1.0 / double( extraNodesPerSegment + 1 );
683
684 int nPoints = ring->numPoints();
685 outX.reserve( ( extraNodesPerSegment + 1 ) * nPoints );
686 outY.reserve( ( extraNodesPerSegment + 1 ) * nPoints );
687 bool withZ = ring->is3D();
688 if ( withZ )
689 outZ.reserve( ( extraNodesPerSegment + 1 ) * nPoints );
690 bool withM = ring->isMeasure();
691 if ( withM )
692 outM.reserve( ( extraNodesPerSegment + 1 ) * nPoints );
693 double x1 = 0;
694 double x2 = 0;
695 double y1 = 0;
696 double y2 = 0;
697 double z1 = 0;
698 double z2 = 0;
699 double m1 = 0;
700 double m2 = 0;
701 double xOut = 0;
702 double yOut = 0;
703 double zOut = 0;
704 double mOut = 0;
705 int extraNodesThisSegment = extraNodesPerSegment;
706 for ( int i = 0; i < nPoints - 1; ++i )
707 {
708 x1 = ring->xAt( i );
709 x2 = ring->xAt( i + 1 );
710 y1 = ring->yAt( i );
711 y2 = ring->yAt( i + 1 );
712 if ( withZ )
713 {
714 z1 = ring->zAt( i );
715 z2 = ring->zAt( i + 1 );
716 }
717 if ( withM )
718 {
719 m1 = ring->mAt( i );
720 m2 = ring->mAt( i + 1 );
721 }
722
723 outX << x1;
724 outY << y1;
725 if ( withZ )
726 outZ << z1;
727 if ( withM )
728 outM << m1;
729
730 if ( extraNodesPerSegment < 0 )
731 {
732 // distance mode
733 extraNodesThisSegment = std::floor( std::sqrt( ( x2 - x1 ) * ( x2 - x1 ) + ( y2 - y1 ) * ( y2 - y1 ) ) / distance );
734 if ( extraNodesThisSegment >= 1 )
735 multiplier = 1.0 / ( extraNodesThisSegment + 1 );
736 }
737
738 for ( int j = 0; j < extraNodesThisSegment; ++j )
739 {
740 double delta = multiplier * ( j + 1 );
741 xOut = x1 + delta * ( x2 - x1 );
742 yOut = y1 + delta * ( y2 - y1 );
743 if ( withZ )
744 zOut = z1 + delta * ( z2 - z1 );
745 if ( withM )
746 mOut = m1 + delta * ( m2 - m1 );
747
748 outX << xOut;
749 outY << yOut;
750 if ( withZ )
751 outZ << zOut;
752 if ( withM )
753 outM << mOut;
754 }
755 }
756 outX << ring->xAt( nPoints - 1 );
757 outY << ring->yAt( nPoints - 1 );
758 if ( withZ )
759 outZ << ring->zAt( nPoints - 1 );
760 if ( withM )
761 outM << ring->mAt( nPoints - 1 );
762
763 QgsLineString *result = new QgsLineString( outX, outY, outZ, outM );
764 return result;
765}
766
767QgsAbstractGeometry *densifyGeometry( const QgsAbstractGeometry *geom, int extraNodesPerSegment = 1, double distance = 1 )
768{
769 std::unique_ptr< QgsAbstractGeometry > segmentizedCopy;
770 if ( QgsWkbTypes::isCurvedType( geom->wkbType() ) )
771 {
772 segmentizedCopy.reset( geom->segmentize() );
773 geom = segmentizedCopy.get();
774 }
775
777 {
778 return doDensify( static_cast< const QgsLineString * >( geom ), extraNodesPerSegment, distance );
779 }
780 else
781 {
782 // polygon
783 const QgsPolygon *polygon = static_cast< const QgsPolygon * >( geom );
784 QgsPolygon *result = new QgsPolygon();
785
786 result->setExteriorRing( doDensify( static_cast< const QgsLineString * >( polygon->exteriorRing() ),
787 extraNodesPerSegment, distance ) );
788 for ( int i = 0; i < polygon->numInteriorRings(); ++i )
789 {
790 result->addInteriorRing( doDensify( static_cast< const QgsLineString * >( polygon->interiorRing( i ) ),
791 extraNodesPerSegment, distance ) );
792 }
793
794 return result;
795 }
796}
797
799{
800 mLastError.clear();
801 if ( !mGeometry )
802 {
803 return QgsGeometry();
804 }
805
807 {
808 return QgsGeometry( mGeometry->clone() ); // point geometry, nothing to do
809 }
810
811 if ( const QgsGeometryCollection *gc = qgsgeometry_cast< const QgsGeometryCollection *>( mGeometry ) )
812 {
813 int numGeom = gc->numGeometries();
814 QVector< QgsAbstractGeometry * > geometryList;
815 geometryList.reserve( numGeom );
816 for ( int i = 0; i < numGeom; ++i )
817 {
818 geometryList << densifyGeometry( gc->geometryN( i ), extraNodesPerSegment );
819 }
820
821 QgsGeometry first = QgsGeometry( geometryList.takeAt( 0 ) );
822 for ( QgsAbstractGeometry *g : std::as_const( geometryList ) )
823 {
824 first.addPart( g );
825 }
826 return first;
827 }
828 else
829 {
830 return QgsGeometry( densifyGeometry( mGeometry, extraNodesPerSegment ) );
831 }
832}
833
835{
836 mLastError.clear();
837 if ( !mGeometry )
838 {
839 return QgsGeometry();
840 }
841
843 {
844 return QgsGeometry( mGeometry->clone() ); // point geometry, nothing to do
845 }
846
847 if ( const QgsGeometryCollection *gc = qgsgeometry_cast< const QgsGeometryCollection *>( mGeometry ) )
848 {
849 int numGeom = gc->numGeometries();
850 QVector< QgsAbstractGeometry * > geometryList;
851 geometryList.reserve( numGeom );
852 for ( int i = 0; i < numGeom; ++i )
853 {
854 geometryList << densifyGeometry( gc->geometryN( i ), -1, distance );
855 }
856
857 QgsGeometry first = QgsGeometry( geometryList.takeAt( 0 ) );
858 for ( QgsAbstractGeometry *g : std::as_const( geometryList ) )
859 {
860 first.addPart( g );
861 }
862 return first;
863 }
864 else
865 {
866 return QgsGeometry( densifyGeometry( mGeometry, -1, distance ) );
867 }
868}
869
871//
872// QgsLineSegmentDistanceComparer
873//
874
875// adapted for QGIS geometry classes from original work at https://github.com/trylock/visibility by trylock
876bool QgsLineSegmentDistanceComparer::operator()( QgsLineSegment2D ab, QgsLineSegment2D cd ) const
877{
878 Q_ASSERT_X( ab.pointLeftOfLine( mOrigin ) != 0,
879 "line_segment_dist_comparer",
880 "AB must not be collinear with the origin." );
881 Q_ASSERT_X( cd.pointLeftOfLine( mOrigin ) != 0,
882 "line_segment_dist_comparer",
883 "CD must not be collinear with the origin." );
884
885 // flip the segments so that if there are common endpoints,
886 // they will be the segment's start points
887 // cppcheck-suppress mismatchingContainerExpression
888 if ( ab.end() == cd.start() || ab.end() == cd.end() )
889 ab.reverse();
890 // cppcheck-suppress mismatchingContainerExpression
891 if ( ab.start() == cd.end() )
892 cd.reverse();
893
894 // cases with common endpoints
895 if ( ab.start() == cd.start() )
896 {
897 const int oad = QgsGeometryUtils::leftOfLine( cd.endX(), cd.endY(), mOrigin.x(), mOrigin.y(), ab.startX(), ab.startY() );
898 const int oab = ab.pointLeftOfLine( mOrigin );
899 // cppcheck-suppress mismatchingContainerExpression
900 if ( ab.end() == cd.end() || oad != oab )
901 return false;
902 else
903 return ab.pointLeftOfLine( cd.end() ) != oab;
904 }
905 else
906 {
907 // cases without common endpoints
908 const int cda = cd.pointLeftOfLine( ab.start() );
909 const int cdb = cd.pointLeftOfLine( ab.end() );
910 if ( cdb == 0 && cda == 0 )
911 {
912 return mOrigin.sqrDist( ab.start() ) < mOrigin.sqrDist( cd.start() );
913 }
914 else if ( cda == cdb || cda == 0 || cdb == 0 )
915 {
916 const int cdo = cd.pointLeftOfLine( mOrigin );
917 return cdo == cda || cdo == cdb;
918 }
919 else
920 {
921 const int abo = ab.pointLeftOfLine( mOrigin );
922 return abo != ab.pointLeftOfLine( cd.start() );
923 }
924 }
925}
926
927//
928// QgsClockwiseAngleComparer
929//
930
931bool QgsClockwiseAngleComparer::operator()( const QgsPointXY &a, const QgsPointXY &b ) const
932{
933 const bool aIsLeft = a.x() < mVertex.x();
934 const bool bIsLeft = b.x() < mVertex.x();
935 if ( aIsLeft != bIsLeft )
936 return bIsLeft;
937
938 if ( qgsDoubleNear( a.x(), mVertex.x() ) && qgsDoubleNear( b.x(), mVertex.x() ) )
939 {
940 if ( a.y() >= mVertex.y() || b.y() >= mVertex.y() )
941 {
942 return b.y() < a.y();
943 }
944 else
945 {
946 return a.y() < b.y();
947 }
948 }
949 else
950 {
951 const QgsVector oa = a - mVertex;
952 const QgsVector ob = b - mVertex;
953 const double det = oa.crossProduct( ob );
954 if ( qgsDoubleNear( det, 0.0 ) )
955 {
956 return oa.lengthSquared() < ob.lengthSquared();
957 }
958 else
959 {
960 return det < 0;
961 }
962 }
963}
964
966
967//
968// QgsRay2D
969//
970
971bool QgsRay2D::intersects( const QgsLineSegment2D &segment, QgsPointXY &intersectPoint ) const
972{
973 const QgsVector ao = origin - segment.start();
974 const QgsVector ab = segment.end() - segment.start();
975 const double det = ab.crossProduct( direction );
976 if ( qgsDoubleNear( det, 0.0 ) )
977 {
978 const int abo = segment.pointLeftOfLine( origin );
979 if ( abo != 0 )
980 {
981 return false;
982 }
983 else
984 {
985 const double distA = ao * direction;
986 const double distB = ( origin - segment.end() ) * direction;
987
988 if ( distA > 0 && distB > 0 )
989 {
990 return false;
991 }
992 else
993 {
994 if ( ( distA > 0 ) != ( distB > 0 ) )
995 intersectPoint = origin;
996 else if ( distA > distB ) // at this point, both distances are negative
997 intersectPoint = segment.start(); // hence the nearest point is A
998 else
999 intersectPoint = segment.end();
1000 return true;
1001 }
1002 }
1003 }
1004 else
1005 {
1006 const double u = ao.crossProduct( direction ) / det;
1007 if ( u < 0.0 || 1.0 < u )
1008 {
1009 return false;
1010 }
1011 else
1012 {
1013 const double t = -ab.crossProduct( ao ) / det;
1014 intersectPoint = origin + direction * t;
1015 return qgsDoubleNear( t, 0.0 ) || t > 0;
1016 }
1017 }
1018}
1019
1020QVector<QgsPointXY> generateSegmentCurve( const QgsPoint &center1, const double radius1, const QgsPoint &center2, const double radius2 )
1021{
1022 // ensure that first circle is smaller than second
1023 if ( radius1 > radius2 )
1024 return generateSegmentCurve( center2, radius2, center1, radius1 );
1025
1026 QgsPointXY t1;
1027 QgsPointXY t2;
1028 QgsPointXY t3;
1029 QgsPointXY t4;
1030 QVector<QgsPointXY> points;
1031 if ( QgsGeometryUtils::circleCircleOuterTangents( center1, radius1, center2, radius2, t1, t2, t3, t4 ) )
1032 {
1033 points << t1
1034 << t2
1035 << t4
1036 << t3;
1037 }
1038 return points;
1039}
1040
1041// partially ported from JTS VariableWidthBuffer,
1042// https://github.com/topobyte/jts/blob/master/jts-lab/src/main/java/com/vividsolutions/jts/operation/buffer/VariableWidthBuffer.java
1043
1044QgsGeometry QgsInternalGeometryEngine::variableWidthBuffer( int segments, const std::function< std::unique_ptr< double[] >( const QgsLineString *line ) > &widthFunction ) const
1045{
1046 mLastError.clear();
1047 if ( !mGeometry )
1048 {
1049 return QgsGeometry();
1050 }
1051
1052 std::vector< std::unique_ptr<QgsLineString > > linesToProcess;
1053
1054 const QgsMultiCurve *multiCurve = qgsgeometry_cast< const QgsMultiCurve * >( mGeometry );
1055 if ( multiCurve )
1056 {
1057 for ( int i = 0; i < multiCurve->partCount(); ++i )
1058 {
1059 if ( static_cast< const QgsCurve * >( multiCurve->geometryN( i ) )->nCoordinates() == 0 )
1060 continue; // skip 0 length lines
1061
1062 linesToProcess.emplace_back( static_cast<QgsLineString *>( multiCurve->geometryN( i )->clone() ) );
1063 }
1064 }
1065
1066 const QgsCurve *curve = qgsgeometry_cast< const QgsCurve * >( mGeometry );
1067 if ( curve )
1068 {
1069 if ( curve->nCoordinates() > 0 )
1070 linesToProcess.emplace_back( static_cast<QgsLineString *>( curve->segmentize() ) );
1071 }
1072
1073 if ( linesToProcess.empty() )
1074 {
1075 QgsGeometry g;
1076 g.mLastError = QStringLiteral( "Input geometry was not a curve type geometry" );
1077 return g;
1078 }
1079
1080 QVector<QgsGeometry> bufferedLines;
1081 bufferedLines.reserve( linesToProcess.size() );
1082
1083 for ( std::unique_ptr< QgsLineString > &line : linesToProcess )
1084 {
1085 QVector<QgsGeometry> parts;
1086 QgsPoint prevPoint;
1087 double prevRadius = 0;
1088 QgsGeometry prevCircle;
1089
1090 std::unique_ptr< double[] > widths = widthFunction( line.get() ) ;
1091 for ( int i = 0; i < line->nCoordinates(); ++i )
1092 {
1093 QgsPoint thisPoint = line->pointN( i );
1094 QgsGeometry thisCircle;
1095 double thisRadius = widths[ i ] / 2.0;
1096 if ( thisRadius > 0 )
1097 {
1098 QgsGeometry p = QgsGeometry( thisPoint.clone() );
1099
1100 QgsCircle circ( thisPoint, thisRadius );
1101 thisCircle = QgsGeometry( circ.toPolygon( segments * 4 ) );
1102 parts << thisCircle;
1103 }
1104 else
1105 {
1106 thisCircle = QgsGeometry( thisPoint.clone() );
1107 }
1108
1109 if ( i > 0 )
1110 {
1111 if ( prevRadius > 0 || thisRadius > 0 )
1112 {
1113 QVector< QgsPointXY > points = generateSegmentCurve( prevPoint, prevRadius, thisPoint, thisRadius );
1114 if ( !points.empty() )
1115 {
1116 // snap points to circle vertices
1117
1118 int atVertex = 0;
1119 int beforeVertex = 0;
1120 int afterVertex = 0;
1121 double sqrDist = 0;
1122 double sqrDistPrev = 0;
1123 for ( int j = 0; j < points.count(); ++j )
1124 {
1125 QgsPointXY pA = prevCircle.closestVertex( points.at( j ), atVertex, beforeVertex, afterVertex, sqrDistPrev );
1126 QgsPointXY pB = thisCircle.closestVertex( points.at( j ), atVertex, beforeVertex, afterVertex, sqrDist );
1127 points[j] = sqrDistPrev < sqrDist ? pA : pB;
1128 }
1129 // close ring
1130 points.append( points.at( 0 ) );
1131
1132 std::unique_ptr< QgsPolygon > poly = std::make_unique< QgsPolygon >();
1133 poly->setExteriorRing( new QgsLineString( points ) );
1134 if ( poly->area() > 0 )
1135 parts << QgsGeometry( std::move( poly ) );
1136 }
1137 }
1138 }
1139 prevPoint = thisPoint;
1140 prevRadius = thisRadius;
1141 prevCircle = thisCircle;
1142 }
1143
1144 bufferedLines << QgsGeometry::unaryUnion( parts );
1145 }
1146
1147 QgsGeometry res = QgsGeometry::collectGeometry( bufferedLines );
1148 // happens on some GEOS versions...
1150 return res;
1151}
1152
1153QgsGeometry QgsInternalGeometryEngine::taperedBuffer( double start, double end, int segments ) const
1154{
1155 mLastError.clear();
1156 start = std::fabs( start );
1157 end = std::fabs( end );
1158
1159 auto interpolateWidths = [ start, end ]( const QgsLineString * line )->std::unique_ptr< double [] >
1160 {
1161 // ported from JTS VariableWidthBuffer,
1162 // https://github.com/topobyte/jts/blob/master/jts-lab/src/main/java/com/vividsolutions/jts/operation/buffer/VariableWidthBuffer.java
1163 std::unique_ptr< double [] > widths( new double[ line->nCoordinates() ] );
1164 widths[0] = start;
1165 widths[line->nCoordinates() - 1] = end;
1166
1167 double lineLength = line->length();
1168 double currentLength = 0;
1169 QgsPoint prevPoint = line->pointN( 0 );
1170 for ( int i = 1; i < line->nCoordinates() - 1; ++i )
1171 {
1172 QgsPoint point = line->pointN( i );
1173 double segmentLength = point.distance( prevPoint );
1174 currentLength += segmentLength;
1175 double lengthFraction = lineLength > 0 ? currentLength / lineLength : 1;
1176 double delta = lengthFraction * ( end - start );
1177 widths[i] = start + delta;
1178 prevPoint = point;
1179 }
1180 return widths;
1181 };
1182
1183 return variableWidthBuffer( segments, interpolateWidths );
1184}
1185
1187{
1188 mLastError.clear();
1189 auto widthByM = []( const QgsLineString * line )->std::unique_ptr< double [] >
1190 {
1191 std::unique_ptr< double [] > widths( new double[ line->nCoordinates() ] );
1192 for ( int i = 0; i < line->nCoordinates(); ++i )
1193 {
1194 widths[ i ] = line->mAt( i );
1195 }
1196 return widths;
1197 };
1198
1199 return variableWidthBuffer( segments, widthByM );
1200}
1201
1202QVector<QgsPointXY> QgsInternalGeometryEngine::randomPointsInPolygon( const QgsGeometry &polygon, int count,
1203 const std::function< bool( const QgsPointXY & ) > &acceptPoint, unsigned long seed, QgsFeedback *feedback, int maxTriesPerPoint )
1204{
1205 if ( polygon.type() != Qgis::GeometryType::Polygon || count == 0 )
1206 return QVector< QgsPointXY >();
1207
1208 // step 1 - tessellate the polygon to triangles
1209 QgsRectangle bounds = polygon.boundingBox();
1210 QgsTessellator t( bounds, false, false, false, true );
1211
1212 if ( polygon.isMultipart() )
1213 {
1214 const QgsMultiSurface *ms = qgsgeometry_cast< const QgsMultiSurface * >( polygon.constGet() );
1215 for ( int i = 0; i < ms->numGeometries(); ++i )
1216 {
1217 if ( feedback && feedback->isCanceled() )
1218 return QVector< QgsPointXY >();
1219
1220 if ( QgsPolygon *poly = qgsgeometry_cast< QgsPolygon * >( ms->geometryN( i ) ) )
1221 {
1222 t.addPolygon( *poly, 0 );
1223 }
1224 else
1225 {
1226 std::unique_ptr< QgsPolygon > p( qgsgeometry_cast< QgsPolygon * >( ms->geometryN( i )->segmentize() ) );
1227 t.addPolygon( *p, 0 );
1228 }
1229 }
1230 }
1231 else
1232 {
1233 if ( const QgsPolygon *poly = qgsgeometry_cast< const QgsPolygon * >( polygon.constGet() ) )
1234 {
1235 t.addPolygon( *poly, 0 );
1236 }
1237 else
1238 {
1239 std::unique_ptr< QgsPolygon > p( qgsgeometry_cast< QgsPolygon * >( polygon.constGet()->segmentize() ) );
1240 t.addPolygon( *p, 0 );
1241 }
1242 }
1243
1244 if ( feedback && feedback->isCanceled() )
1245 return QVector< QgsPointXY >();
1246
1247 const QVector<float> triangleData = t.data();
1248 if ( triangleData.empty() )
1249 return QVector< QgsPointXY >(); //hm
1250
1251 // calculate running sum of triangle areas
1252 std::vector< double > cumulativeAreas;
1253 cumulativeAreas.reserve( t.dataVerticesCount() / 3 );
1254 double totalArea = 0;
1255 for ( auto it = triangleData.constBegin(); it != triangleData.constEnd(); )
1256 {
1257 if ( feedback && feedback->isCanceled() )
1258 return QVector< QgsPointXY >();
1259
1260 const float aX = *it++;
1261 ( void )it++; // z
1262 const float aY = -( *it++ );
1263 const float bX = *it++;
1264 ( void )it++; // z
1265 const float bY = -( *it++ );
1266 const float cX = *it++;
1267 ( void )it++; // z
1268 const float cY = -( *it++ );
1269
1270 const double area = QgsGeometryUtils::triangleArea( aX, aY, bX, bY, cX, cY );
1271 totalArea += area;
1272 cumulativeAreas.emplace_back( totalArea );
1273 }
1274
1275 std::random_device rd;
1276 std::mt19937 mt( seed == 0 ? rd() : seed );
1277 std::uniform_real_distribution<> uniformDist( 0, 1 );
1278
1279 // selects a random triangle, weighted by triangle area
1280 auto selectRandomTriangle = [&cumulativeAreas, totalArea]( double random )->int
1281 {
1282 int triangle = 0;
1283 const double target = random * totalArea;
1284 for ( auto it = cumulativeAreas.begin(); it != cumulativeAreas.end(); ++it, triangle++ )
1285 {
1286 if ( *it > target )
1287 return triangle;
1288 }
1289 Q_ASSERT_X( false, "QgsInternalGeometryEngine::randomPointsInPolygon", "Invalid random triangle index" );
1290 return 0; // no warnings
1291 };
1292
1293
1294 QVector<QgsPointXY> result;
1295 result.reserve( count );
1296 int tries = 0;
1297 for ( int i = 0; i < count; )
1298 {
1299 if ( feedback && feedback->isCanceled() )
1300 return QVector< QgsPointXY >();
1301
1302 const double triangleIndexRnd = uniformDist( mt );
1303 // pick random triangle, weighted by triangle area
1304 const int triangleIndex = selectRandomTriangle( triangleIndexRnd );
1305
1306 // generate a random point inside this triangle
1307 const double weightB = uniformDist( mt );
1308 const double weightC = uniformDist( mt );
1309 double x;
1310 double y;
1311
1312 // get triangle
1313 const double aX = triangleData.at( triangleIndex * 9 ) + bounds.xMinimum();
1314 const double aY = -triangleData.at( triangleIndex * 9 + 2 ) + bounds.yMinimum();
1315 const double bX = triangleData.at( triangleIndex * 9 + 3 ) + bounds.xMinimum();
1316 const double bY = -triangleData.at( triangleIndex * 9 + 5 ) + bounds.yMinimum();
1317 const double cX = triangleData.at( triangleIndex * 9 + 6 ) + bounds.xMinimum();
1318 const double cY = -triangleData.at( triangleIndex * 9 + 8 ) + bounds.yMinimum();
1319 QgsGeometryUtils::weightedPointInTriangle( aX, aY, bX, bY, cX, cY, weightB, weightC, x, y );
1320
1321 QgsPointXY candidate( x, y );
1322 if ( acceptPoint( candidate ) )
1323 {
1324 result << QgsPointXY( x, y );
1325 i++;
1326 tries = 0;
1327 }
1328 else if ( maxTriesPerPoint != 0 )
1329 {
1330 tries++;
1331 // Skip this point if maximum tries is reached
1332 if ( tries == maxTriesPerPoint )
1333 {
1334 tries = 0;
1335 i++;
1336 }
1337 }
1338 }
1339 return result;
1340}
1341
1342// ported from PostGIS' lwgeom pta_unstroke
1343
1344std::unique_ptr< QgsCurve > lineToCurve( const QgsCurve *curve, double distanceTolerance,
1345 double pointSpacingAngleTolerance )
1346{
1347 const Qgis::WkbType flatType = QgsWkbTypes::flatType( curve->wkbType() );
1348 if ( flatType == Qgis::WkbType::CircularString )
1349 {
1350 // already curved, so just return a copy
1351 std::unique_ptr< QgsCircularString > out;
1352 out.reset( qgsgeometry_cast< QgsCircularString * >( curve )->clone() );
1353 return out;
1354 }
1355 else if ( flatType == Qgis::WkbType::CompoundCurve )
1356 {
1357 std::unique_ptr< QgsCompoundCurve > out = std::make_unique< QgsCompoundCurve >();
1358 const QgsCompoundCurve *in = qgsgeometry_cast< const QgsCompoundCurve * >( curve );
1359 for ( int i = 0; i < in->nCurves(); i ++ )
1360 {
1361 std::unique_ptr< QgsCurve > processed = lineToCurve( in->curveAt( i ), distanceTolerance, pointSpacingAngleTolerance );
1362 if ( processed )
1363 {
1364 if ( const QgsCompoundCurve *processedCompoundCurve = qgsgeometry_cast< const QgsCompoundCurve *>( processed.get() ) )
1365 {
1366 for ( int i = 0; i < processedCompoundCurve->nCurves(); ++i )
1367 {
1368 out->addCurve( processedCompoundCurve->curveAt( i )->clone() );
1369 }
1370 }
1371 else
1372 {
1373 out->addCurve( processed.release() );
1374 }
1375 }
1376 }
1377 return out;
1378 }
1379 else if ( flatType == Qgis::WkbType::LineString )
1380 {
1381 const QgsLineString *lineString = qgsgeometry_cast< QgsLineString * >( curve );
1382
1383 std::unique_ptr< QgsCompoundCurve > out = std::make_unique< QgsCompoundCurve >();
1384
1385 /* Minimum number of edges, per quadrant, required to define an arc */
1386 const unsigned int minQuadEdges = 2;
1387
1388 /* Die on null input */
1389 if ( !lineString )
1390 return nullptr;
1391
1392 /* Null on empty input? */
1393 if ( lineString->nCoordinates() == 0 )
1394 return nullptr;
1395
1396 /* We can't desegmentize anything shorter than four points */
1397 if ( lineString->nCoordinates() < 4 )
1398 {
1399 out->addCurve( lineString->clone() );
1400 return out;
1401 }
1402
1403 /* Allocate our result array of vertices that are part of arcs */
1404 int numEdges = lineString->nCoordinates() - 1;
1405 QVector< int > edgesInArcs( numEdges + 1, 0 );
1406
1407 auto arcAngle = []( const QgsPoint & a, const QgsPoint & b, const QgsPoint & c )->double
1408 {
1409 double abX = b.x() - a.x();
1410 double abY = b.y() - a.y();
1411
1412 double cbX = b.x() - c.x();
1413 double cbY = b.y() - c.y();
1414
1415 double dot = ( abX * cbX + abY * cbY ); /* dot product */
1416 double cross = ( abX * cbY - abY * cbX ); /* cross product */
1417
1418 double alpha = std::atan2( cross, dot );
1419
1420 return alpha;
1421 };
1422
1423 /* We make a candidate arc of the first two edges, */
1424 /* And then see if the next edge follows it */
1425 int i = 0;
1426 int j = 0;
1427 int k = 0;
1428 int currentArc = 1;
1429 QgsPoint a1;
1430 QgsPoint a2;
1431 QgsPoint a3;
1432 QgsPoint b;
1433 double centerX = 0.0;
1434 double centerY = 0.0;
1435 double radius = 0;
1436
1437 while ( i < numEdges - 2 )
1438 {
1439 unsigned int arcEdges = 0;
1440 double numQuadrants = 0;
1441 double angle;
1442
1443 bool foundArc = false;
1444 /* Make candidate arc */
1445 a1 = lineString->pointN( i );
1446 a2 = lineString->pointN( i + 1 );
1447 a3 = lineString->pointN( i + 2 );
1448 QgsPoint first = a1;
1449
1450 for ( j = i + 3; j < numEdges + 1; j++ )
1451 {
1452 b = lineString->pointN( j );
1453
1454 /* Does this point fall on our candidate arc? */
1455 if ( QgsGeometryUtils::pointContinuesArc( a1, a2, a3, b, distanceTolerance, pointSpacingAngleTolerance ) )
1456 {
1457 /* Yes. Mark this edge and the two preceding it as arc components */
1458 foundArc = true;
1459 for ( k = j - 1; k > j - 4; k-- )
1460 edgesInArcs[k] = currentArc;
1461 }
1462 else
1463 {
1464 /* No. So we're done with this candidate arc */
1465 currentArc++;
1466 break;
1467 }
1468
1469 a1 = a2;
1470 a2 = a3;
1471 a3 = b;
1472 }
1473 /* Jump past all the edges that were added to the arc */
1474 if ( foundArc )
1475 {
1476 /* Check if an arc was composed by enough edges to be
1477 * really considered an arc
1478 * See http://trac.osgeo.org/postgis/ticket/2420
1479 */
1480 arcEdges = j - 1 - i;
1481 if ( first.x() == b.x() && first.y() == b.y() )
1482 {
1483 numQuadrants = 4;
1484 }
1485 else
1486 {
1487 QgsGeometryUtils::circleCenterRadius( first, b, a1, radius, centerX, centerY );
1488
1489 angle = arcAngle( first, QgsPoint( centerX, centerY ), b );
1490 int p2Side = QgsGeometryUtils::leftOfLine( b.x(), b.y(), first.x(), first.y(), a1.x(), a1.y() );
1491 if ( p2Side >= 0 )
1492 angle = -angle;
1493
1494 if ( angle < 0 )
1495 angle = 2 * M_PI + angle;
1496 numQuadrants = ( 4 * angle ) / ( 2 * M_PI );
1497 }
1498 /* a1 is first point, b is last point */
1499 if ( arcEdges < minQuadEdges * numQuadrants )
1500 {
1501 // LWDEBUGF( 4, "Not enough edges for a %g quadrants arc, %g needed", num_quadrants, min_quad_edges * num_quadrants );
1502 for ( k = j - 1; k >= i; k-- )
1503 edgesInArcs[k] = 0;
1504 }
1505
1506 i = j - 1;
1507 }
1508 else
1509 {
1510 /* Mark this edge as a linear edge */
1511 edgesInArcs[i] = 0;
1512 i = i + 1;
1513 }
1514 }
1515
1516 int start = 0;
1517 int end = 0;
1518 /* non-zero if edge is part of an arc */
1519 int edgeType = edgesInArcs[0];
1520
1521 auto addPointsToCurve = [ lineString, &out ]( int start, int end, int type )
1522 {
1523 if ( type == 0 )
1524 {
1525 // straight segment
1526 QVector< QgsPoint > points;
1527 for ( int j = start; j < end + 2; ++ j )
1528 {
1529 points.append( lineString->pointN( j ) );
1530 }
1531 std::unique_ptr< QgsCurve > straightSegment = std::make_unique< QgsLineString >( points );
1532 out->addCurve( straightSegment.release() );
1533 }
1534 else
1535 {
1536 // curved segment
1537 QVector< QgsPoint > points;
1538 points.append( lineString->pointN( start ) );
1539 points.append( lineString->pointN( ( start + end + 1 ) / 2 ) );
1540 points.append( lineString->pointN( end + 1 ) );
1541 std::unique_ptr< QgsCircularString > curvedSegment = std::make_unique< QgsCircularString >();
1542 curvedSegment->setPoints( points );
1543 out->addCurve( curvedSegment.release() );
1544 }
1545 };
1546
1547 for ( int i = 1; i < numEdges; i++ )
1548 {
1549 if ( edgeType != edgesInArcs[i] )
1550 {
1551 end = i - 1;
1552 addPointsToCurve( start, end, edgeType );
1553 start = i;
1554 edgeType = edgesInArcs[i];
1555 }
1556 }
1557
1558 /* Roll out last item */
1559 end = numEdges - 1;
1560 addPointsToCurve( start, end, edgeType );
1561
1562 // return a simplified type if it doesn't need to be a compound curve resu
1563 if ( QgsWkbTypes::flatType( out->simplifiedTypeRef()->wkbType() ) == Qgis::WkbType::CircularString )
1564 {
1565 std::unique_ptr< QgsCircularString > res;
1566 res.reset( qgsgeometry_cast< const QgsCircularString * >( out->simplifiedTypeRef() )->clone() );
1567 return res;
1568 }
1569
1570 return out;
1571 }
1572
1573 return nullptr;
1574}
1575
1576std::unique_ptr< QgsAbstractGeometry > convertGeometryToCurves( const QgsAbstractGeometry *geom, double distanceTolerance, double angleTolerance )
1577{
1579 {
1580 return lineToCurve( qgsgeometry_cast< const QgsCurve * >( geom ), distanceTolerance, angleTolerance );
1581 }
1582 else
1583 {
1584 // polygon
1585 const QgsCurvePolygon *polygon = qgsgeometry_cast< const QgsCurvePolygon * >( geom );
1586 std::unique_ptr< QgsCurvePolygon > result = std::make_unique< QgsCurvePolygon>();
1587
1588 result->setExteriorRing( lineToCurve( polygon->exteriorRing(),
1589 distanceTolerance, angleTolerance ).release() );
1590 for ( int i = 0; i < polygon->numInteriorRings(); ++i )
1591 {
1592 result->addInteriorRing( lineToCurve( polygon->interiorRing( i ),
1593 distanceTolerance, angleTolerance ).release() );
1594 }
1595
1596 return result;
1597 }
1598}
1599
1600QgsGeometry QgsInternalGeometryEngine::convertToCurves( double distanceTolerance, double angleTolerance ) const
1601{
1602 mLastError.clear();
1603 if ( !mGeometry )
1604 {
1605 return QgsGeometry();
1606 }
1607
1609 {
1610 return QgsGeometry( mGeometry->clone() ); // point geometry, nothing to do
1611 }
1612
1613 if ( const QgsGeometryCollection *gc = qgsgeometry_cast< const QgsGeometryCollection *>( mGeometry ) )
1614 {
1615 int numGeom = gc->numGeometries();
1616 QVector< QgsAbstractGeometry * > geometryList;
1617 geometryList.reserve( numGeom );
1618 for ( int i = 0; i < numGeom; ++i )
1619 {
1620 geometryList << convertGeometryToCurves( gc->geometryN( i ), distanceTolerance, angleTolerance ).release();
1621 }
1622
1623 QgsGeometry first = QgsGeometry( geometryList.takeAt( 0 ) );
1624 for ( QgsAbstractGeometry *g : std::as_const( geometryList ) )
1625 {
1626 first.addPart( g );
1627 }
1628 return first;
1629 }
1630 else
1631 {
1632 return QgsGeometry( convertGeometryToCurves( mGeometry, distanceTolerance, angleTolerance ) );
1633 }
1634}
1635
1636QgsGeometry QgsInternalGeometryEngine::orientedMinimumBoundingBox( double &area, double &angle, double &width, double &height ) const
1637{
1638 mLastError.clear();
1639
1640 QgsRectangle minRect;
1641 area = std::numeric_limits<double>::max();
1642 angle = 0;
1643 width = std::numeric_limits<double>::max();
1644 height = std::numeric_limits<double>::max();
1645
1646 if ( !mGeometry || mGeometry->nCoordinates() < 2 )
1647 return QgsGeometry();
1648
1649 std::unique_ptr< QgsGeometryEngine >engine( QgsGeometry::createGeometryEngine( mGeometry ) );
1650 std::unique_ptr< QgsAbstractGeometry > hull( engine->convexHull( &mLastError ) );
1651 if ( !hull )
1652 return QgsGeometry();
1653
1654 QgsVertexId vertexId;
1655 QgsPoint pt0;
1656 QgsPoint pt1;
1657 QgsPoint pt2;
1658 // get first point
1659 hull->nextVertex( vertexId, pt0 );
1660 pt1 = pt0;
1661 double totalRotation = 0;
1662 while ( hull->nextVertex( vertexId, pt2 ) )
1663 {
1664 double currentAngle = QgsGeometryUtils::lineAngle( pt1.x(), pt1.y(), pt2.x(), pt2.y() );
1665 double rotateAngle = 180.0 / M_PI * currentAngle;
1666 totalRotation += rotateAngle;
1667
1668 QTransform t = QTransform::fromTranslate( pt0.x(), pt0.y() );
1669 t.rotate( rotateAngle );
1670 t.translate( -pt0.x(), -pt0.y() );
1671
1672 hull->transform( t );
1673
1674 QgsRectangle bounds = hull->boundingBox();
1675 double currentArea = bounds.width() * bounds.height();
1676 if ( currentArea < area )
1677 {
1678 minRect = bounds;
1679 area = currentArea;
1680 angle = totalRotation;
1681 width = bounds.width();
1682 height = bounds.height();
1683 }
1684
1685 pt1 = hull->vertexAt( vertexId );
1686 }
1687
1688 QgsGeometry minBounds = QgsGeometry::fromRect( minRect );
1689 minBounds.rotate( angle, QgsPointXY( pt0.x(), pt0.y() ) );
1690
1691 if ( width > height )
1692 {
1693 width = minRect.height();
1694 height = minRect.width();
1695 angle = angle + 90.0;
1696 }
1697
1698 // constrain angle to 0 - 180
1699 if ( angle > 180.0 )
1700 angle = std::fmod( angle, 180.0 );
1701
1702 return minBounds;
1703}
1704
1705std::unique_ptr< QgsLineString > triangularWavesAlongLine( const QgsLineString *line, double wavelength, const double amplitude, const bool strictWavelength )
1706{
1707 const int totalPoints = line->numPoints();
1708 if ( totalPoints < 2 )
1709 return nullptr;
1710
1711 const double *x = line->xData();
1712 const double *y = line->yData();
1713
1714 double prevX = *x++;
1715 double prevY = *y++;
1716
1717 QVector< double > outX;
1718 QVector< double > outY;
1719 const double totalLength = line->length();
1720
1721 const int maxRepetitions = std::ceil( totalLength / wavelength );
1722 if ( !strictWavelength )
1723 wavelength = totalLength / maxRepetitions;
1724
1725 const int estimatedPoints = maxRepetitions * 2 + 2;
1726 outX.reserve( estimatedPoints );
1727 outY.reserve( estimatedPoints );
1728 outX.append( prevX );
1729 outY.append( prevY );
1730
1731 double distanceToNextPointFromStartOfSegment = wavelength / 4;
1732 int side = -1;
1733 for ( int i = 1; i < totalPoints; ++i )
1734 {
1735 double thisX = *x++;
1736 double thisY = *y++;
1737
1738 const double segmentAngleRadians = QgsGeometryUtils::lineAngle( prevX, prevY, thisX, thisY );
1739 const double segmentLength = std::sqrt( ( thisX - prevX ) * ( thisX - prevX ) + ( thisY - prevY ) * ( thisY - prevY ) );
1740 while ( distanceToNextPointFromStartOfSegment < segmentLength || qgsDoubleNear( distanceToNextPointFromStartOfSegment, segmentLength ) )
1741 {
1742 // point falls on this segment - truncate to segment length if qgsDoubleNear test was actually > segment length
1743 const double distanceToPoint = std::min( distanceToNextPointFromStartOfSegment, segmentLength );
1744 double pX, pY;
1745 QgsGeometryUtils::pointOnLineWithDistance( prevX, prevY, thisX, thisY, distanceToPoint, pX, pY );
1746
1747 // project point on line out by amplitude
1748 const double outPointX = pX + side * amplitude * std::sin( segmentAngleRadians + M_PI_2 );
1749 const double outPointY = pY + side * amplitude * std::cos( segmentAngleRadians + M_PI_2 );
1750
1751 outX.append( outPointX );
1752 outY.append( outPointY );
1753
1754 distanceToNextPointFromStartOfSegment += wavelength / 2;
1755 side = -side;
1756 }
1757
1758 prevX = thisX;
1759 prevY = thisY;
1760 distanceToNextPointFromStartOfSegment -= segmentLength;
1761 }
1762
1763 outX.append( prevX );
1764 outY.append( prevY );
1765
1766 return std::make_unique< QgsLineString >( outX, outY );
1767}
1768
1769std::unique_ptr< QgsLineString > triangularWavesRandomizedAlongLine( const QgsLineString *line,
1770 const double minimumWavelength, const double maximumWavelength,
1771 const double minimumAmplitude, const double maximumAmplitude,
1772 std::uniform_real_distribution<> &uniformDist, std::mt19937 &mt )
1773{
1774 const int totalPoints = line->numPoints();
1775 if ( totalPoints < 2 )
1776 return nullptr;
1777
1778 const double *x = line->xData();
1779 const double *y = line->yData();
1780
1781 double prevX = *x++;
1782 double prevY = *y++;
1783
1784 QVector< double > outX;
1785 QVector< double > outY;
1786 const double totalLength = line->length();
1787
1788 const int maxRepetitions = std::ceil( totalLength / minimumWavelength );
1789
1790 const int estimatedPoints = maxRepetitions * 2 + 2;
1791 outX.reserve( estimatedPoints );
1792 outY.reserve( estimatedPoints );
1793 outX.append( prevX );
1794 outY.append( prevY );
1795
1796 double wavelength = uniformDist( mt ) * ( maximumWavelength - minimumWavelength ) + minimumWavelength;
1797 double distanceToNextPointFromStartOfSegment = wavelength / 4;
1798 double amplitude = uniformDist( mt ) * ( maximumAmplitude - minimumAmplitude ) + minimumAmplitude;
1799
1800 int side = -1;
1801 for ( int i = 1; i < totalPoints; ++i )
1802 {
1803 double thisX = *x++;
1804 double thisY = *y++;
1805
1806 const double segmentAngleRadians = QgsGeometryUtils::lineAngle( prevX, prevY, thisX, thisY );
1807 const double segmentLength = std::sqrt( ( thisX - prevX ) * ( thisX - prevX ) + ( thisY - prevY ) * ( thisY - prevY ) );
1808 while ( distanceToNextPointFromStartOfSegment < segmentLength || qgsDoubleNear( distanceToNextPointFromStartOfSegment, segmentLength ) )
1809 {
1810 // point falls on this segment - truncate to segment length if qgsDoubleNear test was actually > segment length
1811 const double distanceToPoint = std::min( distanceToNextPointFromStartOfSegment, segmentLength );
1812 double pX, pY;
1813 QgsGeometryUtils::pointOnLineWithDistance( prevX, prevY, thisX, thisY, distanceToPoint, pX, pY );
1814
1815 // project point on line out by amplitude
1816 const double outPointX = pX + side * amplitude * std::sin( segmentAngleRadians + M_PI_2 );
1817 const double outPointY = pY + side * amplitude * std::cos( segmentAngleRadians + M_PI_2 );
1818 outX.append( outPointX );
1819 outY.append( outPointY );
1820
1821 wavelength = uniformDist( mt ) * ( maximumWavelength - minimumWavelength ) + minimumWavelength;
1822 amplitude = uniformDist( mt ) * ( maximumAmplitude - minimumAmplitude ) + minimumAmplitude;
1823
1824 distanceToNextPointFromStartOfSegment += wavelength / 2;
1825 side = -side;
1826 }
1827
1828 prevX = thisX;
1829 prevY = thisY;
1830 distanceToNextPointFromStartOfSegment -= segmentLength;
1831 }
1832
1833 outX.append( prevX );
1834 outY.append( prevY );
1835
1836 return std::make_unique< QgsLineString >( outX, outY );
1837}
1838
1839std::unique_ptr< QgsAbstractGeometry > triangularWavesPrivate( const QgsAbstractGeometry *geom, double wavelength, double amplitude, bool strictWavelength )
1840{
1841 std::unique_ptr< QgsAbstractGeometry > segmentizedCopy;
1842 if ( QgsWkbTypes::isCurvedType( geom->wkbType() ) )
1843 {
1844 segmentizedCopy.reset( geom->segmentize() );
1845 geom = segmentizedCopy.get();
1846 }
1847
1849 {
1850 return triangularWavesAlongLine( static_cast< const QgsLineString * >( geom ), wavelength, amplitude, strictWavelength );
1851 }
1852 else
1853 {
1854 // polygon
1855 const QgsPolygon *polygon = static_cast< const QgsPolygon * >( geom );
1856 std::unique_ptr< QgsPolygon > result = std::make_unique< QgsPolygon >();
1857
1858 result->setExteriorRing( triangularWavesAlongLine( static_cast< const QgsLineString * >( polygon->exteriorRing() ),
1859 wavelength, amplitude, strictWavelength ).release() );
1860 for ( int i = 0; i < polygon->numInteriorRings(); ++i )
1861 {
1862 result->addInteriorRing( triangularWavesAlongLine( static_cast< const QgsLineString * >( polygon->interiorRing( i ) ),
1863 wavelength, amplitude, strictWavelength ).release() );
1864 }
1865
1866 return result;
1867 }
1868}
1869
1870std::unique_ptr< QgsAbstractGeometry > triangularWavesRandomizedPrivate( const QgsAbstractGeometry *geom, double minimumWavelength, double maximumWavelength, double minimumAmplitude, double maximumAmplitude, std::uniform_real_distribution<> &uniformDist, std::mt19937 &mt )
1871{
1872 std::unique_ptr< QgsAbstractGeometry > segmentizedCopy;
1873 if ( QgsWkbTypes::isCurvedType( geom->wkbType() ) )
1874 {
1875 segmentizedCopy.reset( geom->segmentize() );
1876 geom = segmentizedCopy.get();
1877 }
1878
1880 {
1881 return triangularWavesRandomizedAlongLine( static_cast< const QgsLineString * >( geom ), minimumWavelength, maximumWavelength, minimumAmplitude, maximumAmplitude, uniformDist, mt );
1882 }
1883 else
1884 {
1885 // polygon
1886 const QgsPolygon *polygon = static_cast< const QgsPolygon * >( geom );
1887 std::unique_ptr< QgsPolygon > result = std::make_unique< QgsPolygon >();
1888
1889 result->setExteriorRing( triangularWavesRandomizedAlongLine( static_cast< const QgsLineString * >( polygon->exteriorRing() ),
1890 minimumWavelength, maximumWavelength, minimumAmplitude, maximumAmplitude, uniformDist, mt ).release() );
1891 for ( int i = 0; i < polygon->numInteriorRings(); ++i )
1892 {
1893 result->addInteriorRing( triangularWavesRandomizedAlongLine( static_cast< const QgsLineString * >( polygon->interiorRing( i ) ),
1894 minimumWavelength, maximumWavelength, minimumAmplitude, maximumAmplitude, uniformDist, mt ).release() );
1895 }
1896
1897 return result;
1898 }
1899}
1900
1901QgsGeometry QgsInternalGeometryEngine::triangularWaves( double wavelength, double amplitude, bool strictWavelength ) const
1902{
1903 if ( wavelength < 0 || qgsDoubleNear( wavelength, 0 ) )
1904 return QgsGeometry();
1905
1906 mLastError.clear();
1907 if ( !mGeometry )
1908 {
1909 return QgsGeometry();
1910 }
1911
1913 {
1914 return QgsGeometry( mGeometry->clone() ); // point geometry, nothing to do
1915 }
1916
1917 if ( const QgsGeometryCollection *gc = qgsgeometry_cast< const QgsGeometryCollection *>( mGeometry ) )
1918 {
1919 int numGeom = gc->numGeometries();
1920 QVector< QgsAbstractGeometry * > geometryList;
1921 geometryList.reserve( numGeom );
1922 for ( int i = 0; i < numGeom; ++i )
1923 {
1924 geometryList << triangularWavesPrivate( gc->geometryN( i ), wavelength, amplitude, strictWavelength ).release();
1925 }
1926
1927 QgsGeometry first = QgsGeometry( geometryList.takeAt( 0 ) );
1928 for ( QgsAbstractGeometry *g : std::as_const( geometryList ) )
1929 {
1930 first.addPart( g );
1931 }
1932 return first;
1933 }
1934 else
1935 {
1936 return QgsGeometry( triangularWavesPrivate( mGeometry, wavelength, amplitude, strictWavelength ) );
1937 }
1938}
1939
1940QgsGeometry QgsInternalGeometryEngine::triangularWavesRandomized( double minimumWavelength, double maximumWavelength, double minimumAmplitude, double maximumAmplitude, unsigned long seed ) const
1941{
1942 if ( minimumWavelength < 0 || qgsDoubleNear( minimumWavelength, 0 ) || maximumWavelength < 0 || qgsDoubleNear( maximumWavelength, 0 ) || maximumWavelength < minimumWavelength )
1943 return QgsGeometry();
1944
1945 mLastError.clear();
1946 if ( !mGeometry )
1947 {
1948 return QgsGeometry();
1949 }
1950
1951 std::random_device rd;
1952 std::mt19937 mt( seed == 0 ? rd() : seed );
1953 std::uniform_real_distribution<> uniformDist( 0, 1 );
1954
1956 {
1957 return QgsGeometry( mGeometry->clone() ); // point geometry, nothing to do
1958 }
1959
1960 if ( const QgsGeometryCollection *gc = qgsgeometry_cast< const QgsGeometryCollection *>( mGeometry ) )
1961 {
1962 int numGeom = gc->numGeometries();
1963 QVector< QgsAbstractGeometry * > geometryList;
1964 geometryList.reserve( numGeom );
1965 for ( int i = 0; i < numGeom; ++i )
1966 {
1967 geometryList << triangularWavesRandomizedPrivate( gc->geometryN( i ), minimumWavelength, maximumWavelength, minimumAmplitude, maximumAmplitude, uniformDist, mt ).release();
1968 }
1969
1970 QgsGeometry first = QgsGeometry( geometryList.takeAt( 0 ) );
1971 for ( QgsAbstractGeometry *g : std::as_const( geometryList ) )
1972 {
1973 first.addPart( g );
1974 }
1975 return first;
1976 }
1977 else
1978 {
1979 return QgsGeometry( triangularWavesRandomizedPrivate( mGeometry, minimumWavelength, maximumWavelength, minimumAmplitude, maximumAmplitude, uniformDist, mt ) );
1980 }
1981}
1982
1983std::unique_ptr< QgsLineString > squareWavesAlongLine( const QgsLineString *line, double wavelength, const double amplitude, const bool strictWavelength )
1984{
1985 const int totalPoints = line->numPoints();
1986 if ( totalPoints < 2 )
1987 return nullptr;
1988
1989 const double *x = line->xData();
1990 const double *y = line->yData();
1991
1992 double prevX = *x++;
1993 double prevY = *y++;
1994
1995 QVector< double > outX;
1996 QVector< double > outY;
1997 const double totalLength = line->length();
1998
1999 const int maxRepetitions = std::ceil( totalLength / wavelength );
2000 if ( !strictWavelength )
2001 wavelength = totalLength / maxRepetitions;
2002
2003 const int estimatedPoints = maxRepetitions * 4 + 2;
2004 outX.reserve( estimatedPoints );
2005 outY.reserve( estimatedPoints );
2006 outX.append( prevX );
2007 outY.append( prevY );
2008
2009 const double startSegmentAngleRadians = QgsGeometryUtils::lineAngle( prevX, prevY, *x, *y );
2010 outX.append( prevX - amplitude * std::sin( startSegmentAngleRadians + M_PI_2 ) );
2011 outY.append( prevY - amplitude * std::cos( startSegmentAngleRadians + M_PI_2 ) );
2012
2013 double distanceToNextPointFromStartOfSegment = wavelength / 2;
2014
2015 int side = -1;
2016 for ( int i = 1; i < totalPoints; ++i )
2017 {
2018 double thisX = *x++;
2019 double thisY = *y++;
2020
2021 const double segmentAngleRadians = QgsGeometryUtils::lineAngle( prevX, prevY, thisX, thisY );
2022 const double segmentLength = std::sqrt( ( thisX - prevX ) * ( thisX - prevX ) + ( thisY - prevY ) * ( thisY - prevY ) );
2023 while ( distanceToNextPointFromStartOfSegment < segmentLength || qgsDoubleNear( distanceToNextPointFromStartOfSegment, segmentLength ) )
2024 {
2025 // point falls on this segment - truncate to segment length if qgsDoubleNear test was actually > segment length
2026 const double distanceToPoint = std::min( distanceToNextPointFromStartOfSegment, segmentLength );
2027 double pX, pY;
2028 QgsGeometryUtils::pointOnLineWithDistance( prevX, prevY, thisX, thisY, distanceToPoint, pX, pY );
2029
2030 // project point on line out by amplitude
2031 const double sinAngle = std::sin( segmentAngleRadians + M_PI_2 );
2032 const double cosAngle = std::cos( segmentAngleRadians + M_PI_2 );
2033 outX.append( pX + side * amplitude * sinAngle );
2034 outY.append( pY + side * amplitude * cosAngle );
2035 outX.append( pX - side * amplitude * sinAngle );
2036 outY.append( pY - side * amplitude * cosAngle );
2037
2038 distanceToNextPointFromStartOfSegment += wavelength / 2;
2039 side = -side;
2040 }
2041
2042 prevX = thisX;
2043 prevY = thisY;
2044 distanceToNextPointFromStartOfSegment -= segmentLength;
2045 }
2046
2047 // replace last point, which will be projected out to amplitude of wave, with the actual end point of the line
2048 outX.pop_back();
2049 outY.pop_back();
2050 outX.append( prevX );
2051 outY.append( prevY );
2052
2053 return std::make_unique< QgsLineString >( outX, outY );
2054}
2055
2056std::unique_ptr< QgsLineString > squareWavesRandomizedAlongLine( const QgsLineString *line,
2057 const double minimumWavelength, const double maximumWavelength,
2058 const double minimumAmplitude, const double maximumAmplitude,
2059 std::uniform_real_distribution<> &uniformDist, std::mt19937 &mt )
2060{
2061 const int totalPoints = line->numPoints();
2062 if ( totalPoints < 2 )
2063 return nullptr;
2064
2065 const double *x = line->xData();
2066 const double *y = line->yData();
2067
2068 double prevX = *x++;
2069 double prevY = *y++;
2070
2071 QVector< double > outX;
2072 QVector< double > outY;
2073 const double totalLength = line->length();
2074
2075 const int maxRepetitions = std::ceil( totalLength / minimumWavelength );
2076
2077 const int estimatedPoints = maxRepetitions * 4 + 2;
2078 outX.reserve( estimatedPoints );
2079 outY.reserve( estimatedPoints );
2080 outX.append( prevX );
2081 outY.append( prevY );
2082
2083 double amplitude = uniformDist( mt ) * ( maximumAmplitude - minimumAmplitude ) + minimumAmplitude;
2084
2085 double segmentAngleRadians = QgsGeometryUtils::lineAngle( prevX, prevY, *x, *y );
2086 outX.append( prevX - amplitude * std::sin( segmentAngleRadians + M_PI_2 ) );
2087 outY.append( prevY - amplitude * std::cos( segmentAngleRadians + M_PI_2 ) );
2088
2089 double wavelength = uniformDist( mt ) * ( maximumWavelength - minimumWavelength ) + minimumWavelength;
2090 double distanceToNextPointFromStartOfSegment = wavelength / 2;
2091
2092 int side = -1;
2093 for ( int i = 1; i < totalPoints; ++i )
2094 {
2095 double thisX = *x++;
2096 double thisY = *y++;
2097
2098 segmentAngleRadians = QgsGeometryUtils::lineAngle( prevX, prevY, thisX, thisY );
2099 const double segmentLength = std::sqrt( ( thisX - prevX ) * ( thisX - prevX ) + ( thisY - prevY ) * ( thisY - prevY ) );
2100 while ( distanceToNextPointFromStartOfSegment < segmentLength || qgsDoubleNear( distanceToNextPointFromStartOfSegment, segmentLength ) )
2101 {
2102 // point falls on this segment - truncate to segment length if qgsDoubleNear test was actually > segment length
2103 const double distanceToPoint = std::min( distanceToNextPointFromStartOfSegment, segmentLength );
2104 double pX, pY;
2105 QgsGeometryUtils::pointOnLineWithDistance( prevX, prevY, thisX, thisY, distanceToPoint, pX, pY );
2106
2107 // project point on line out by amplitude
2108 const double sinAngle = std::sin( segmentAngleRadians + M_PI_2 );
2109 const double cosAngle = std::cos( segmentAngleRadians + M_PI_2 );
2110 outX.append( pX + side * amplitude * sinAngle );
2111 outY.append( pY + side * amplitude * cosAngle );
2112
2113 amplitude = uniformDist( mt ) * ( maximumAmplitude - minimumAmplitude ) + minimumAmplitude;
2114 outX.append( pX - side * amplitude * sinAngle );
2115 outY.append( pY - side * amplitude * cosAngle );
2116
2117 wavelength = uniformDist( mt ) * ( maximumWavelength - minimumWavelength ) + minimumWavelength;
2118 distanceToNextPointFromStartOfSegment += wavelength / 2;
2119 side = -side;
2120 }
2121
2122 prevX = thisX;
2123 prevY = thisY;
2124 distanceToNextPointFromStartOfSegment -= segmentLength;
2125 }
2126
2127 outX.append( prevX + side * amplitude * std::sin( segmentAngleRadians + M_PI_2 ) );
2128 outY.append( prevY + side * amplitude * std::cos( segmentAngleRadians + M_PI_2 ) );
2129 outX.append( prevX );
2130 outY.append( prevY );
2131
2132 return std::make_unique< QgsLineString >( outX, outY );
2133}
2134
2135std::unique_ptr< QgsAbstractGeometry > squareWavesPrivate( const QgsAbstractGeometry *geom, double wavelength, double amplitude, bool strictWavelength )
2136{
2137 std::unique_ptr< QgsAbstractGeometry > segmentizedCopy;
2138 if ( QgsWkbTypes::isCurvedType( geom->wkbType() ) )
2139 {
2140 segmentizedCopy.reset( geom->segmentize() );
2141 geom = segmentizedCopy.get();
2142 }
2143
2145 {
2146 return squareWavesAlongLine( static_cast< const QgsLineString * >( geom ), wavelength, amplitude, strictWavelength );
2147 }
2148 else
2149 {
2150 // polygon
2151 const QgsPolygon *polygon = static_cast< const QgsPolygon * >( geom );
2152 std::unique_ptr< QgsPolygon > result = std::make_unique< QgsPolygon >();
2153
2154 result->setExteriorRing( squareWavesAlongLine( static_cast< const QgsLineString * >( polygon->exteriorRing() ),
2155 wavelength, amplitude, strictWavelength ).release() );
2156 for ( int i = 0; i < polygon->numInteriorRings(); ++i )
2157 {
2158 result->addInteriorRing( squareWavesAlongLine( static_cast< const QgsLineString * >( polygon->interiorRing( i ) ),
2159 wavelength, amplitude, strictWavelength ).release() );
2160 }
2161
2162 return result;
2163 }
2164}
2165
2166std::unique_ptr< QgsAbstractGeometry > squareWavesRandomizedPrivate( const QgsAbstractGeometry *geom, double minimumWavelength, double maximumWavelength, double minimumAmplitude, double maximumAmplitude, std::uniform_real_distribution<> &uniformDist, std::mt19937 &mt )
2167{
2168 std::unique_ptr< QgsAbstractGeometry > segmentizedCopy;
2169 if ( QgsWkbTypes::isCurvedType( geom->wkbType() ) )
2170 {
2171 segmentizedCopy.reset( geom->segmentize() );
2172 geom = segmentizedCopy.get();
2173 }
2174
2176 {
2177 return squareWavesRandomizedAlongLine( static_cast< const QgsLineString * >( geom ), minimumWavelength, maximumWavelength, minimumAmplitude, maximumAmplitude, uniformDist, mt );
2178 }
2179 else
2180 {
2181 // polygon
2182 const QgsPolygon *polygon = static_cast< const QgsPolygon * >( geom );
2183 std::unique_ptr< QgsPolygon > result = std::make_unique< QgsPolygon >();
2184
2185 result->setExteriorRing( squareWavesRandomizedAlongLine( static_cast< const QgsLineString * >( polygon->exteriorRing() ),
2186 minimumWavelength, maximumWavelength, minimumAmplitude, maximumAmplitude, uniformDist, mt ).release() );
2187 for ( int i = 0; i < polygon->numInteriorRings(); ++i )
2188 {
2189 result->addInteriorRing( squareWavesRandomizedAlongLine( static_cast< const QgsLineString * >( polygon->interiorRing( i ) ),
2190 minimumWavelength, maximumWavelength, minimumAmplitude, maximumAmplitude, uniformDist, mt ).release() );
2191 }
2192
2193 return result;
2194 }
2195}
2196
2197QgsGeometry QgsInternalGeometryEngine::squareWaves( double wavelength, double amplitude, bool strictWavelength ) const
2198{
2199 if ( wavelength < 0 || qgsDoubleNear( wavelength, 0 ) )
2200 return QgsGeometry();
2201
2202 mLastError.clear();
2203 if ( !mGeometry )
2204 {
2205 return QgsGeometry();
2206 }
2207
2209 {
2210 return QgsGeometry( mGeometry->clone() ); // point geometry, nothing to do
2211 }
2212
2213 if ( const QgsGeometryCollection *gc = qgsgeometry_cast< const QgsGeometryCollection *>( mGeometry ) )
2214 {
2215 int numGeom = gc->numGeometries();
2216 QVector< QgsAbstractGeometry * > geometryList;
2217 geometryList.reserve( numGeom );
2218 for ( int i = 0; i < numGeom; ++i )
2219 {
2220 geometryList << squareWavesPrivate( gc->geometryN( i ), wavelength, amplitude, strictWavelength ).release();
2221 }
2222
2223 QgsGeometry first = QgsGeometry( geometryList.takeAt( 0 ) );
2224 for ( QgsAbstractGeometry *g : std::as_const( geometryList ) )
2225 {
2226 first.addPart( g );
2227 }
2228 return first;
2229 }
2230 else
2231 {
2232 return QgsGeometry( squareWavesPrivate( mGeometry, wavelength, amplitude, strictWavelength ) );
2233 }
2234}
2235
2236QgsGeometry QgsInternalGeometryEngine::squareWavesRandomized( double minimumWavelength, double maximumWavelength, double minimumAmplitude, double maximumAmplitude, unsigned long seed ) const
2237{
2238 if ( minimumWavelength < 0 || qgsDoubleNear( minimumWavelength, 0 ) || maximumWavelength < 0 || qgsDoubleNear( maximumWavelength, 0 ) || maximumWavelength < minimumWavelength )
2239 return QgsGeometry();
2240
2241 mLastError.clear();
2242 if ( !mGeometry )
2243 {
2244 return QgsGeometry();
2245 }
2246
2247 std::random_device rd;
2248 std::mt19937 mt( seed == 0 ? rd() : seed );
2249 std::uniform_real_distribution<> uniformDist( 0, 1 );
2250
2252 {
2253 return QgsGeometry( mGeometry->clone() ); // point geometry, nothing to do
2254 }
2255
2256 if ( const QgsGeometryCollection *gc = qgsgeometry_cast< const QgsGeometryCollection *>( mGeometry ) )
2257 {
2258 int numGeom = gc->numGeometries();
2259 QVector< QgsAbstractGeometry * > geometryList;
2260 geometryList.reserve( numGeom );
2261 for ( int i = 0; i < numGeom; ++i )
2262 {
2263 geometryList << squareWavesRandomizedPrivate( gc->geometryN( i ), minimumWavelength, maximumWavelength, minimumAmplitude, maximumAmplitude, uniformDist, mt ).release();
2264 }
2265
2266 QgsGeometry first = QgsGeometry( geometryList.takeAt( 0 ) );
2267 for ( QgsAbstractGeometry *g : std::as_const( geometryList ) )
2268 {
2269 first.addPart( g );
2270 }
2271 return first;
2272 }
2273 else
2274 {
2275 return QgsGeometry( squareWavesRandomizedPrivate( mGeometry, minimumWavelength, maximumWavelength, minimumAmplitude, maximumAmplitude, uniformDist, mt ) );
2276 }
2277}
2278
2279std::unique_ptr< QgsLineString > roundWavesAlongLine( const QgsLineString *line, double wavelength, const double amplitude, const bool strictWavelength )
2280{
2281 const int totalPoints = line->numPoints();
2282 if ( totalPoints < 2 )
2283 return nullptr;
2284
2285 const double *x = line->xData();
2286 const double *y = line->yData();
2287
2288 double prevX = *x++;
2289 double prevY = *y++;
2290
2291 const double totalLength = line->length();
2292
2293 const int maxRepetitions = std::ceil( totalLength / wavelength );
2294 if ( !strictWavelength )
2295 wavelength = totalLength / maxRepetitions;
2296
2297 const int segments = 10;
2298
2299 int side = -1;
2300
2301 double xOutBuffer[4] { prevX, prevX, prevX, prevX };
2302 double yOutBuffer[4] { prevY, prevY, prevY, prevY };
2303 bool isFirstPart = true;
2304
2305 double distanceToNextPointFromStartOfSegment = wavelength / 8;
2306 int bufferIndex = 1;
2307 std::unique_ptr< QgsLineString > out = std::make_unique< QgsLineString >();
2308
2309 double segmentAngleRadians = 0;
2310 double remainingDistance = totalLength;
2311 double totalCoveredDistance = 0;
2312
2313 for ( int i = 1; i < totalPoints; ++i )
2314 {
2315 double thisX = *x++;
2316 double thisY = *y++;
2317
2318 segmentAngleRadians = QgsGeometryUtils::lineAngle( prevX, prevY, thisX, thisY );
2319 const double segmentLength = std::sqrt( ( thisX - prevX ) * ( thisX - prevX ) + ( thisY - prevY ) * ( thisY - prevY ) );
2320 while ( distanceToNextPointFromStartOfSegment < segmentLength || qgsDoubleNear( distanceToNextPointFromStartOfSegment, segmentLength ) )
2321 {
2322 // point falls on this segment - truncate to segment length if qgsDoubleNear test was actually > segment length
2323 const double distanceToPoint = std::min( distanceToNextPointFromStartOfSegment, segmentLength );
2324 double pX, pY;
2325 QgsGeometryUtils::pointOnLineWithDistance( prevX, prevY, thisX, thisY, distanceToPoint, pX, pY );
2326 remainingDistance = totalLength - totalCoveredDistance - distanceToPoint;
2327
2328 const double sinAngle = std::sin( segmentAngleRadians + M_PI_2 );
2329 const double cosAngle = std::cos( segmentAngleRadians + M_PI_2 );
2330
2331 if ( bufferIndex == 0 )
2332 {
2333 xOutBuffer[0] = pX + side * amplitude * sinAngle;
2334 yOutBuffer[0] = pY + side * amplitude * cosAngle;
2335 bufferIndex = 1;
2336 distanceToNextPointFromStartOfSegment += wavelength / 4;
2337 }
2338 else if ( bufferIndex == 1 && isFirstPart )
2339 {
2340 xOutBuffer[1] = ( xOutBuffer[0] + pX - side * amplitude * sinAngle ) * 0.5;
2341 yOutBuffer[1] = ( yOutBuffer[0] + pY - side * amplitude * cosAngle ) * 0.5;
2342 xOutBuffer[2] = pX - side * amplitude * sinAngle;
2343 yOutBuffer[2] = pY - side * amplitude * cosAngle;
2344 bufferIndex = 2;
2345 distanceToNextPointFromStartOfSegment += wavelength / 8;
2346 }
2347 else if ( bufferIndex == 1 )
2348 {
2349 xOutBuffer[1] = pX + side * amplitude * sinAngle;
2350 yOutBuffer[1] = pY + side * amplitude * cosAngle;
2351 xOutBuffer[2] = pX - side * amplitude * sinAngle;
2352 yOutBuffer[2] = pY - side * amplitude * cosAngle;
2353 bufferIndex = 2;
2354 distanceToNextPointFromStartOfSegment += wavelength / 4;
2355 }
2356 else if ( bufferIndex == 2 )
2357 {
2358 xOutBuffer[3] = pX - side * amplitude * sinAngle;
2359 yOutBuffer[3] = pY - side * amplitude * cosAngle;
2360
2361 if ( isFirstPart )
2362 {
2363 xOutBuffer[2] = ( xOutBuffer[2] + xOutBuffer[3] ) * 0.5;
2364 yOutBuffer[2] = ( yOutBuffer[2] + yOutBuffer[3] ) * 0.5;
2365 isFirstPart = false;
2366 }
2367
2368 // flush curve
2369 std::unique_ptr< QgsLineString > bezier( QgsLineString::fromBezierCurve( QgsPoint( xOutBuffer[0], yOutBuffer[0] ),
2370 QgsPoint( xOutBuffer[1], yOutBuffer[1] ),
2371 QgsPoint( xOutBuffer[2], yOutBuffer[2] ),
2372 QgsPoint( xOutBuffer[3], yOutBuffer[3] ),
2373 segments ) );
2374 out->append( bezier.get() );
2375
2376 // shuffle buffer alone
2377 xOutBuffer[0] = xOutBuffer[3];
2378 yOutBuffer[0] = yOutBuffer[3];
2379 bufferIndex = 1;
2380 side = -side;
2381 distanceToNextPointFromStartOfSegment += wavelength / 4;
2382 }
2383 }
2384 totalCoveredDistance += segmentLength;
2385 prevX = thisX;
2386 prevY = thisY;
2387 distanceToNextPointFromStartOfSegment -= segmentLength;
2388 }
2389
2390 const double midX = prevX - remainingDistance / 2 * std::sin( segmentAngleRadians );
2391 const double midY = prevY - remainingDistance / 2 * std::cos( segmentAngleRadians );
2392 const double pX = midX + side * amplitude * std::sin( segmentAngleRadians + M_PI_2 );
2393 const double pY = midY + side * amplitude * std::cos( segmentAngleRadians + M_PI_2 );
2394
2395 std::unique_ptr< QgsLineString > bezier( QgsLineString::fromBezierCurve( out->endPoint(),
2396 QgsPoint( ( out->endPoint().x() + pX ) * 0.5, ( out->endPoint().y() + pY ) * 0.5 ),
2397 QgsPoint( ( prevX + pX ) * 0.5, ( prevY + pY ) * 0.5 ),
2398 QgsPoint( prevX, prevY ),
2399 segments / 2 ) );
2400 out->append( bezier.get() );
2401
2402 return out;
2403}
2404
2405std::unique_ptr< QgsLineString > roundWavesRandomizedAlongLine( const QgsLineString *line,
2406 const double minimumWavelength, const double maximumWavelength,
2407 const double minimumAmplitude, const double maximumAmplitude,
2408 std::uniform_real_distribution<> &uniformDist, std::mt19937 &mt )
2409{
2410 const int totalPoints = line->numPoints();
2411 if ( totalPoints < 2 )
2412 return nullptr;
2413
2414 const double *x = line->xData();
2415 const double *y = line->yData();
2416
2417 double prevX = *x++;
2418 double prevY = *y++;
2419
2420 const double totalLength = line->length();
2421
2422 const int segments = 10;
2423
2424 int side = -1;
2425
2426 double xOutBuffer[4] { prevX, prevX, prevX, prevX };
2427 double yOutBuffer[4] { prevY, prevY, prevY, prevY };
2428 bool isFirstPart = true;
2429
2430 double amplitude = uniformDist( mt ) * ( maximumAmplitude - minimumAmplitude ) + minimumAmplitude;
2431 double wavelength = uniformDist( mt ) * ( maximumWavelength - minimumWavelength ) + minimumWavelength;
2432
2433 double distanceToNextPointFromStartOfSegment = wavelength / 8;
2434 int bufferIndex = 1;
2435 std::unique_ptr< QgsLineString > out = std::make_unique< QgsLineString >();
2436
2437 double segmentAngleRadians = 0;
2438
2439 double remainingDistance = totalLength;
2440 double totalCoveredDistance = 0;
2441
2442 for ( int i = 1; i < totalPoints; ++i )
2443 {
2444 double thisX = *x++;
2445 double thisY = *y++;
2446
2447 segmentAngleRadians = QgsGeometryUtils::lineAngle( prevX, prevY, thisX, thisY );
2448 const double segmentLength = std::sqrt( ( thisX - prevX ) * ( thisX - prevX ) + ( thisY - prevY ) * ( thisY - prevY ) );
2449 while ( distanceToNextPointFromStartOfSegment < segmentLength || qgsDoubleNear( distanceToNextPointFromStartOfSegment, segmentLength ) )
2450 {
2451 // point falls on this segment - truncate to segment length if qgsDoubleNear test was actually > segment length
2452 const double distanceToPoint = std::min( distanceToNextPointFromStartOfSegment, segmentLength );
2453 double pX, pY;
2454 QgsGeometryUtils::pointOnLineWithDistance( prevX, prevY, thisX, thisY, distanceToPoint, pX, pY );
2455 remainingDistance = totalLength - totalCoveredDistance - distanceToPoint;
2456
2457 const double sinAngle = std::sin( segmentAngleRadians + M_PI_2 );
2458 const double cosAngle = std::cos( segmentAngleRadians + M_PI_2 );
2459
2460 if ( bufferIndex == 0 )
2461 {
2462 xOutBuffer[0] = pX + side * amplitude * sinAngle;
2463 yOutBuffer[0] = pY + side * amplitude * cosAngle;
2464 bufferIndex = 1;
2465 distanceToNextPointFromStartOfSegment += wavelength / 4;
2466 }
2467 else if ( bufferIndex == 1 && isFirstPart )
2468 {
2469 xOutBuffer[1] = ( xOutBuffer[0] + pX - side * amplitude * sinAngle ) * 0.5;
2470 yOutBuffer[1] = ( yOutBuffer[0] + pY - side * amplitude * cosAngle ) * 0.5;
2471 xOutBuffer[2] = pX - side * amplitude * sinAngle;
2472 yOutBuffer[2] = pY - side * amplitude * cosAngle;
2473 bufferIndex = 2;
2474 distanceToNextPointFromStartOfSegment += wavelength / 8;
2475 }
2476 else if ( bufferIndex == 1 )
2477 {
2478 xOutBuffer[1] = pX + side * amplitude * sinAngle;
2479 yOutBuffer[1] = pY + side * amplitude * cosAngle;
2480 amplitude = uniformDist( mt ) * ( maximumAmplitude - minimumAmplitude ) + minimumAmplitude;
2481 xOutBuffer[2] = pX - side * amplitude * sinAngle;
2482 bufferIndex = 2;
2483 yOutBuffer[2] = pY - side * amplitude * cosAngle;
2484 distanceToNextPointFromStartOfSegment += wavelength / 4;
2485 }
2486 else if ( bufferIndex == 2 )
2487 {
2488 xOutBuffer[3] = pX - side * amplitude * sinAngle;
2489 yOutBuffer[3] = pY - side * amplitude * cosAngle;
2490
2491 if ( isFirstPart )
2492 {
2493 xOutBuffer[2] = ( xOutBuffer[2] + xOutBuffer[3] ) * 0.5;
2494 yOutBuffer[2] = ( yOutBuffer[2] + yOutBuffer[3] ) * 0.5;
2495 isFirstPart = false;
2496 }
2497
2498 // flush curve
2499 std::unique_ptr< QgsLineString > bezier( QgsLineString::fromBezierCurve( QgsPoint( xOutBuffer[0], yOutBuffer[0] ),
2500 QgsPoint( xOutBuffer[1], yOutBuffer[1] ),
2501 QgsPoint( xOutBuffer[2], yOutBuffer[2] ),
2502 QgsPoint( xOutBuffer[3], yOutBuffer[3] ),
2503 segments ) );
2504 out->append( bezier.get() );
2505
2506 // shuffle buffer alone
2507 xOutBuffer[0] = xOutBuffer[3];
2508 yOutBuffer[0] = yOutBuffer[3];
2509 bufferIndex = 1;
2510 side = -side;
2511
2512 wavelength = uniformDist( mt ) * ( maximumWavelength - minimumWavelength ) + minimumWavelength;
2513
2514 distanceToNextPointFromStartOfSegment += wavelength / 4;
2515 }
2516 }
2517 totalCoveredDistance += segmentLength;
2518
2519 prevX = thisX;
2520 prevY = thisY;
2521 distanceToNextPointFromStartOfSegment -= segmentLength;
2522 }
2523
2524 const double midX = prevX - remainingDistance / 2 * std::sin( segmentAngleRadians );
2525 const double midY = prevY - remainingDistance / 2 * std::cos( segmentAngleRadians );
2526 const double pX = midX + side * amplitude * std::sin( segmentAngleRadians + M_PI_2 );
2527 const double pY = midY + side * amplitude * std::cos( segmentAngleRadians + M_PI_2 );
2528
2529 if ( out->isEmpty() )
2530 {
2531 std::unique_ptr< QgsLineString > bezier( QgsLineString::fromBezierCurve( line->startPoint(),
2532 QgsPoint( ( line->startPoint().x() + pX ) * 0.5, ( line->startPoint().y() + pY ) * 0.5 ),
2533 QgsPoint( ( prevX + pX ) * 0.5, ( prevY + pY ) * 0.5 ),
2534 QgsPoint( prevX, prevY ),
2535 segments ) );
2536 out->append( bezier.get() );
2537 }
2538 else
2539 {
2540 std::unique_ptr< QgsLineString > bezier( QgsLineString::fromBezierCurve( out->endPoint(),
2541 QgsPoint( ( out->endPoint().x() + pX ) * 0.5, ( out->endPoint().y() + pY ) * 0.5 ),
2542 QgsPoint( ( prevX + pX ) * 0.5, ( prevY + pY ) * 0.5 ),
2543 QgsPoint( prevX, prevY ),
2544 segments ) );
2545 out->append( bezier.get() );
2546 }
2547
2548 return out;
2549}
2550
2551std::unique_ptr< QgsAbstractGeometry > roundWavesPrivate( const QgsAbstractGeometry *geom, double wavelength, double amplitude, bool strictWavelength )
2552{
2553 std::unique_ptr< QgsAbstractGeometry > segmentizedCopy;
2554 if ( QgsWkbTypes::isCurvedType( geom->wkbType() ) )
2555 {
2556 segmentizedCopy.reset( geom->segmentize() );
2557 geom = segmentizedCopy.get();
2558 }
2559
2561 {
2562 return roundWavesAlongLine( static_cast< const QgsLineString * >( geom ), wavelength, amplitude, strictWavelength );
2563 }
2564 else
2565 {
2566 // polygon
2567 const QgsPolygon *polygon = static_cast< const QgsPolygon * >( geom );
2568 std::unique_ptr< QgsPolygon > result = std::make_unique< QgsPolygon >();
2569
2570 result->setExteriorRing( roundWavesAlongLine( static_cast< const QgsLineString * >( polygon->exteriorRing() ),
2571 wavelength, amplitude, strictWavelength ).release() );
2572 for ( int i = 0; i < polygon->numInteriorRings(); ++i )
2573 {
2574 result->addInteriorRing( roundWavesAlongLine( static_cast< const QgsLineString * >( polygon->interiorRing( i ) ),
2575 wavelength, amplitude, strictWavelength ).release() );
2576 }
2577
2578 return result;
2579 }
2580}
2581
2582std::unique_ptr< QgsAbstractGeometry > roundWavesRandomizedPrivate( const QgsAbstractGeometry *geom, double minimumWavelength, double maximumWavelength, double minimumAmplitude, double maximumAmplitude, std::uniform_real_distribution<> &uniformDist, std::mt19937 &mt )
2583{
2584 std::unique_ptr< QgsAbstractGeometry > segmentizedCopy;
2585 if ( QgsWkbTypes::isCurvedType( geom->wkbType() ) )
2586 {
2587 segmentizedCopy.reset( geom->segmentize() );
2588 geom = segmentizedCopy.get();
2589 }
2590
2592 {
2593 return roundWavesRandomizedAlongLine( static_cast< const QgsLineString * >( geom ), minimumWavelength, maximumWavelength, minimumAmplitude, maximumAmplitude, uniformDist, mt );
2594 }
2595 else
2596 {
2597 // polygon
2598 const QgsPolygon *polygon = static_cast< const QgsPolygon * >( geom );
2599 std::unique_ptr< QgsPolygon > result = std::make_unique< QgsPolygon >();
2600
2601 result->setExteriorRing( roundWavesRandomizedAlongLine( static_cast< const QgsLineString * >( polygon->exteriorRing() ),
2602 minimumWavelength, maximumWavelength, minimumAmplitude, maximumAmplitude, uniformDist, mt ).release() );
2603 for ( int i = 0; i < polygon->numInteriorRings(); ++i )
2604 {
2605 result->addInteriorRing( roundWavesRandomizedAlongLine( static_cast< const QgsLineString * >( polygon->interiorRing( i ) ),
2606 minimumWavelength, maximumWavelength, minimumAmplitude, maximumAmplitude, uniformDist, mt ).release() );
2607 }
2608
2609 return result;
2610 }
2611}
2612
2613QgsGeometry QgsInternalGeometryEngine::roundWaves( double wavelength, double amplitude, bool strictWavelength ) const
2614{
2615 if ( wavelength < 0 || qgsDoubleNear( wavelength, 0 ) )
2616 return QgsGeometry();
2617
2618 mLastError.clear();
2619 if ( !mGeometry )
2620 {
2621 return QgsGeometry();
2622 }
2623
2625 {
2626 return QgsGeometry( mGeometry->clone() ); // point geometry, nothing to do
2627 }
2628
2629 if ( const QgsGeometryCollection *gc = qgsgeometry_cast< const QgsGeometryCollection *>( mGeometry ) )
2630 {
2631 int numGeom = gc->numGeometries();
2632 QVector< QgsAbstractGeometry * > geometryList;
2633 geometryList.reserve( numGeom );
2634 for ( int i = 0; i < numGeom; ++i )
2635 {
2636 geometryList << roundWavesPrivate( gc->geometryN( i ), wavelength, amplitude, strictWavelength ).release();
2637 }
2638
2639 QgsGeometry first = QgsGeometry( geometryList.takeAt( 0 ) );
2640 for ( QgsAbstractGeometry *g : std::as_const( geometryList ) )
2641 {
2642 first.addPart( g );
2643 }
2644 return first;
2645 }
2646 else
2647 {
2648 return QgsGeometry( roundWavesPrivate( mGeometry, wavelength, amplitude, strictWavelength ) );
2649 }
2650}
2651
2652QgsGeometry QgsInternalGeometryEngine::roundWavesRandomized( double minimumWavelength, double maximumWavelength, double minimumAmplitude, double maximumAmplitude, unsigned long seed ) const
2653{
2654 if ( minimumWavelength < 0 || qgsDoubleNear( minimumWavelength, 0 ) || maximumWavelength < 0 || qgsDoubleNear( maximumWavelength, 0 ) || maximumWavelength < minimumWavelength )
2655 return QgsGeometry();
2656
2657 mLastError.clear();
2658 if ( !mGeometry )
2659 {
2660 return QgsGeometry();
2661 }
2662
2663 std::random_device rd;
2664 std::mt19937 mt( seed == 0 ? rd() : seed );
2665 std::uniform_real_distribution<> uniformDist( 0, 1 );
2666
2668 {
2669 return QgsGeometry( mGeometry->clone() ); // point geometry, nothing to do
2670 }
2671
2672 if ( const QgsGeometryCollection *gc = qgsgeometry_cast< const QgsGeometryCollection *>( mGeometry ) )
2673 {
2674 int numGeom = gc->numGeometries();
2675 QVector< QgsAbstractGeometry * > geometryList;
2676 geometryList.reserve( numGeom );
2677 for ( int i = 0; i < numGeom; ++i )
2678 {
2679 geometryList << roundWavesRandomizedPrivate( gc->geometryN( i ), minimumWavelength, maximumWavelength, minimumAmplitude, maximumAmplitude, uniformDist, mt ).release();
2680 }
2681
2682 QgsGeometry first = QgsGeometry( geometryList.takeAt( 0 ) );
2683 for ( QgsAbstractGeometry *g : std::as_const( geometryList ) )
2684 {
2685 first.addPart( g );
2686 }
2687 return first;
2688 }
2689 else
2690 {
2691 return QgsGeometry( roundWavesRandomizedPrivate( mGeometry, minimumWavelength, maximumWavelength, minimumAmplitude, maximumAmplitude, uniformDist, mt ) );
2692 }
2693}
2694
2695std::unique_ptr< QgsMultiLineString > dashPatternAlongLine( const QgsLineString *line,
2696 const QVector< double> &pattern,
2700 double patternOffset )
2701{
2702 const int totalPoints = line->numPoints();
2703 if ( totalPoints < 2 )
2704 return nullptr;
2705
2706 const int patternSize = pattern.size();
2707
2708 const double *x = line->xData();
2709 const double *y = line->yData();
2710
2711 double prevX = *x++;
2712 double prevY = *y++;
2713
2714 std::unique_ptr< QgsMultiLineString > result = std::make_unique< QgsMultiLineString >();
2715
2716 QVector< double > outX;
2717 QVector< double > outY;
2718 const double totalLength = line->length();
2719
2720 double patternLength = 0;
2721 double patternDashLength = 0;
2722 double patternGapLength = 0;
2723 for ( int i = 0; i < pattern.size(); ++i )
2724 {
2725 patternLength += pattern.at( i );
2726 if ( i % 2 == 0 )
2727 patternDashLength += pattern.at( i );
2728 else
2729 patternGapLength += pattern.at( i );
2730 }
2731
2732 double firstPatternLength = 0;
2733 double firstPatternDashLength = 0;
2734 double firstPatternGapLength = 0;
2735 switch ( startRule )
2736 {
2739 firstPatternLength = patternLength;
2740 firstPatternDashLength = patternDashLength;
2741 firstPatternGapLength = patternGapLength;
2742 break;
2744 firstPatternLength = patternLength - pattern.at( 0 ) * 0.5; // remove half of first dash
2745 firstPatternDashLength = patternDashLength - pattern.at( 0 ) * 0.5;
2746 firstPatternGapLength = patternGapLength;
2747 break;
2749 firstPatternLength = pattern.at( patternSize - 1 ); // start at LAST dash
2750 firstPatternDashLength = 0;
2751 firstPatternGapLength = pattern.at( patternSize - 1 );
2752 break;
2754 firstPatternLength = pattern.at( patternSize - 1 ) * 0.5; // start at half of last dash
2755 firstPatternDashLength = 0;
2756 firstPatternGapLength = pattern.at( patternSize - 1 ) * 0.5;
2757 break;
2758 }
2759
2760 const bool isSmallEnoughForSinglePattern = ( totalLength - firstPatternLength ) < patternLength * 0.5;
2761
2762 double lastPatternLength = isSmallEnoughForSinglePattern ? firstPatternLength : patternLength;
2763 double lastPatternDashLength = isSmallEnoughForSinglePattern ? firstPatternDashLength : patternDashLength;
2764 double lastPatternGapLength = isSmallEnoughForSinglePattern ? firstPatternGapLength : patternGapLength;
2765 switch ( endRule )
2766 {
2768 lastPatternLength = 0;
2769 lastPatternDashLength = 0;
2770 lastPatternGapLength = 0;
2771 break;
2773 lastPatternLength -= pattern.at( patternSize - 1 ); // remove last gap
2774 lastPatternGapLength -= pattern.at( patternSize - 1 );
2775 break;
2777 lastPatternLength -= pattern.at( patternSize - 1 ) + pattern.at( patternSize - 2 ) * 0.5; // remove last gap, half of last dash
2778 lastPatternDashLength -= pattern.at( patternSize - 2 ) * 0.5;
2779 lastPatternGapLength -= pattern.at( patternSize - 1 );
2780 break;
2782 lastPatternGapLength = patternGapLength;
2783 break;
2785 lastPatternLength -= pattern.at( patternSize - 1 ) * 0.5; // remove half of last gap
2786 lastPatternGapLength -= pattern.at( patternSize - 1 ) * 0.5;
2787 break;
2788 }
2789
2790 const double remainingLengthForCompletePatterns = totalLength - ( !isSmallEnoughForSinglePattern ? firstPatternLength : 0 ) - lastPatternLength;
2791 const int middlePatternRepetitions = std::max( static_cast< int >( std::round( remainingLengthForCompletePatterns / patternLength ) ), 0 );
2792
2793 const double totalUnscaledLengthOfPatterns = ( !isSmallEnoughForSinglePattern ? firstPatternLength : 0 ) + middlePatternRepetitions * patternLength + lastPatternLength;
2794 const double totalUnscaledLengthOfDashes = ( !isSmallEnoughForSinglePattern ? firstPatternDashLength : 0 ) + middlePatternRepetitions * patternDashLength + lastPatternDashLength;
2795 const double totalUnscaledLengthOfGaps = ( !isSmallEnoughForSinglePattern ? firstPatternGapLength : 0 ) + middlePatternRepetitions * patternGapLength + lastPatternGapLength;
2796
2797 double dashLengthScalingFactor = 1;
2798 double gapLengthScalingFactor = 1;
2800 {
2801 // calculate scaling factors
2802 const double lengthToShrinkBy = totalUnscaledLengthOfPatterns - totalLength;
2803
2804 switch ( adjustment )
2805 {
2807 dashLengthScalingFactor = totalLength / totalUnscaledLengthOfPatterns;
2808 gapLengthScalingFactor = dashLengthScalingFactor;
2809 break;
2811 dashLengthScalingFactor = ( totalUnscaledLengthOfDashes - lengthToShrinkBy ) / totalUnscaledLengthOfDashes;
2812 break;
2814 gapLengthScalingFactor = ( totalUnscaledLengthOfGaps - lengthToShrinkBy ) / totalUnscaledLengthOfGaps;
2815 break;
2816 }
2817 }
2818
2819 dashLengthScalingFactor = std::max( dashLengthScalingFactor, 0.0 );
2820 gapLengthScalingFactor = std::max( gapLengthScalingFactor, 0.0 );
2821
2822 const int maxPatterns = middlePatternRepetitions + 2;
2823 result->reserve( maxPatterns );
2824
2825 int patternIndex = 0;
2826 double distanceToNextPointFromStartOfSegment = pattern.at( 0 ) * dashLengthScalingFactor;
2827 bool isDash = true;
2828 switch ( startRule )
2829 {
2832 break;
2834 distanceToNextPointFromStartOfSegment *= 0.5; // skip to half way through first dash
2835 break;
2837 patternIndex = patternSize - 1; // start at last gap
2838 isDash = false;
2839 distanceToNextPointFromStartOfSegment = pattern.at( patternSize - 1 ) * gapLengthScalingFactor;
2840 break;
2842 patternIndex = patternSize - 1; // skip straight to half way through last gap
2843 isDash = false;
2844 distanceToNextPointFromStartOfSegment = 0.5 * pattern.at( patternSize - 1 ) * gapLengthScalingFactor;
2845 break;
2846 }
2847
2848 const double adjustedOffset = fmod( patternOffset, patternLength );
2849 const double scaledOffset = ( adjustedOffset < 0 ? ( adjustedOffset + patternLength ) : adjustedOffset ) * ( gapLengthScalingFactor + dashLengthScalingFactor ) / 2;
2850 if ( !qgsDoubleNear( scaledOffset, 0 ) )
2851 {
2852 // shuffle pattern along by offset
2853 double remainingOffset = scaledOffset;
2854 while ( remainingOffset > 0 )
2855 {
2856 if ( distanceToNextPointFromStartOfSegment > remainingOffset )
2857 {
2858 distanceToNextPointFromStartOfSegment -= remainingOffset;
2859 break;
2860 }
2861
2862 remainingOffset -= distanceToNextPointFromStartOfSegment;
2863 isDash = !isDash;
2864 patternIndex++;
2865 if ( patternIndex == patternSize )
2866 patternIndex = 0;
2867
2868 distanceToNextPointFromStartOfSegment = pattern.at( patternIndex ) * ( isDash ? dashLengthScalingFactor : gapLengthScalingFactor );
2869 }
2870 }
2871
2872 if ( isDash )
2873 {
2874 outX.append( prevX );
2875 outY.append( prevY );
2876 }
2877
2878 for ( int i = 1; i < totalPoints; ++i )
2879 {
2880 double thisX = *x++;
2881 double thisY = *y++;
2882
2883 const double segmentLength = std::sqrt( ( thisX - prevX ) * ( thisX - prevX ) + ( thisY - prevY ) * ( thisY - prevY ) );
2884 while ( distanceToNextPointFromStartOfSegment < segmentLength || qgsDoubleNear( distanceToNextPointFromStartOfSegment, segmentLength ) )
2885 {
2886 // point falls on this segment - truncate to segment length if qgsDoubleNear test was actually > segment length
2887 const double distanceToPoint = std::min( distanceToNextPointFromStartOfSegment, segmentLength );
2888 double pX, pY;
2889 QgsGeometryUtils::pointOnLineWithDistance( prevX, prevY, thisX, thisY, distanceToPoint, pX, pY );
2890
2891 outX.append( pX );
2892 outY.append( pY );
2893 if ( isDash )
2894 {
2895 result->addGeometry( new QgsLineString( outX, outY ) );
2896 outX.resize( 0 );
2897 outY.resize( 0 );
2898 }
2899
2900 isDash = !isDash;
2901 patternIndex++;
2902 if ( patternIndex >= patternSize )
2903 patternIndex = 0;
2904
2905 distanceToNextPointFromStartOfSegment += pattern.at( patternIndex ) * ( isDash ? dashLengthScalingFactor : gapLengthScalingFactor );
2906 }
2907
2908 if ( isDash )
2909 {
2910 outX.append( thisX );
2911 outY.append( thisY );
2912 }
2913
2914 prevX = thisX;
2915 prevY = thisY;
2916 distanceToNextPointFromStartOfSegment -= segmentLength;
2917 }
2918
2919 if ( isDash )
2920 {
2921 outX.append( prevX );
2922 outY.append( prevY );
2923 result->addGeometry( new QgsLineString( outX, outY ) );
2924 }
2925
2926 return result;
2927}
2928
2929std::unique_ptr< QgsAbstractGeometry > applyDashPatternPrivate( const QgsAbstractGeometry *geom,
2930 const QVector<double> &pattern,
2934 double patternOffset )
2935{
2936 std::unique_ptr< QgsAbstractGeometry > segmentizedCopy;
2937 if ( QgsWkbTypes::isCurvedType( geom->wkbType() ) )
2938 {
2939 segmentizedCopy.reset( geom->segmentize() );
2940 geom = segmentizedCopy.get();
2941 }
2942
2944 {
2945 return dashPatternAlongLine( static_cast< const QgsLineString * >( geom ), pattern, startRule, endRule, adjustment, patternOffset );
2946 }
2947 else
2948 {
2949 // polygon
2950 const QgsPolygon *polygon = static_cast< const QgsPolygon * >( geom );
2951 std::unique_ptr< QgsMultiLineString > result = std::make_unique< QgsMultiLineString >();
2952
2953 std::unique_ptr< QgsMultiLineString > exteriorParts = dashPatternAlongLine( static_cast< const QgsLineString * >( polygon->exteriorRing() ), pattern, startRule, endRule, adjustment, patternOffset );
2954 for ( int i = 0; i < exteriorParts->numGeometries(); ++i )
2955 result->addGeometry( exteriorParts->geometryN( i )->clone() );
2956
2957 for ( int i = 0; i < polygon->numInteriorRings(); ++i )
2958 {
2959 std::unique_ptr< QgsMultiLineString > ringParts = dashPatternAlongLine( static_cast< const QgsLineString * >( polygon->interiorRing( i ) ), pattern, startRule, endRule, adjustment, patternOffset );
2960 for ( int j = 0; j < ringParts->numGeometries(); ++j )
2961 result->addGeometry( ringParts->geometryN( j )->clone() );
2962 }
2963
2964 return result;
2965 }
2966}
2967
2969{
2970 if ( pattern.size() < 2 )
2971 return QgsGeometry( mGeometry->clone() );
2972
2973 mLastError.clear();
2974 if ( !mGeometry || mGeometry->isEmpty() )
2975 {
2976 return QgsGeometry();
2977 }
2978
2980 {
2981 return QgsGeometry( mGeometry->clone() ); // point geometry, nothing to do
2982 }
2983
2984 if ( const QgsGeometryCollection *gc = qgsgeometry_cast< const QgsGeometryCollection *>( mGeometry ) )
2985 {
2986 int numGeom = gc->numGeometries();
2987 QVector< QgsAbstractGeometry * > geometryList;
2988 geometryList.reserve( numGeom );
2989 for ( int i = 0; i < numGeom; ++i )
2990 {
2991 geometryList << applyDashPatternPrivate( gc->geometryN( i ), pattern, startRule, endRule, adjustment, patternOffset ).release();
2992 }
2993
2994 QgsGeometry first = QgsGeometry( geometryList.takeAt( 0 ) );
2995 for ( QgsAbstractGeometry *g : std::as_const( geometryList ) )
2996 {
2997 if ( const QgsGeometryCollection *collection = qgsgeometry_cast< QgsGeometryCollection * >( g ) )
2998 {
2999 for ( int j = 0; j < collection->numGeometries(); ++j )
3000 {
3001 first.addPart( collection->geometryN( j )->clone() );
3002 }
3003 delete collection;
3004 }
3005 else
3006 {
3007 first.addPart( g );
3008 }
3009 }
3010 return first;
3011 }
3012 else
3013 {
3014 return QgsGeometry( applyDashPatternPrivate( mGeometry, pattern, startRule, endRule, adjustment, patternOffset ) );
3015 }
3016}
DashPatternSizeAdjustment
Dash pattern size adjustment options.
Definition qgis.h:2421
@ ScaleDashOnly
Only dash lengths are adjusted.
@ ScaleBothDashAndGap
Both the dash and gap lengths are adjusted equally.
@ ScaleGapOnly
Only gap lengths are adjusted.
@ Polygon
Polygons.
DashPatternLineEndingRule
Dash pattern line ending rules.
Definition qgis.h:2406
@ HalfDash
Start or finish the pattern with a half length dash.
@ HalfGap
Start or finish the pattern with a half length gap.
@ FullGap
Start or finish the pattern with a full gap.
@ FullDash
Start or finish the pattern with a full dash.
WkbType
The WKB type describes the number of dimensions a geometry has.
Definition qgis.h:182
@ CompoundCurve
CompoundCurve.
@ LineString
LineString.
@ Polygon
Polygon.
@ CircularString
CircularString.
The vertex_iterator class provides STL-style iterator for vertices.
Abstract base class for all geometries.
vertex_iterator vertices_end() const
Returns STL-style iterator pointing to the imaginary vertex after the last vertex of the geometry.
virtual QgsAbstractGeometry * segmentize(double tolerance=M_PI/180., SegmentationToleranceType toleranceType=MaximumAngle) const
Returns a version of the geometry without curves.
bool isMeasure() const
Returns true if the geometry contains m values.
virtual QgsRectangle boundingBox() const
Returns the minimal bounding box for the geometry.
bool is3D() const
Returns true if the geometry is 3D and contains a z-value.
virtual int nCoordinates() const
Returns the number of nodes contained in the geometry.
Qgis::WkbType wkbType() const
Returns the WKB type of the geometry.
virtual bool isEmpty() const
Returns true if the geometry is empty.
vertex_iterator vertices_begin() const
Returns STL-style iterator pointing to the first vertex of the geometry.
virtual QgsAbstractGeometry * clone() const =0
Clones the geometry by performing a deep copy.
Circle geometry type.
Definition qgscircle.h:44
Compound curve geometry type.
int nCurves() const
Returns the number of curves in the geometry.
const QgsCurve * curveAt(int i) const
Returns the curve at the specified index.
Curve polygon geometry type.
int numInteriorRings() const
Returns the number of interior rings contained with the curve polygon.
const QgsCurve * exteriorRing() const
Returns the curve polygon's exterior ring.
const QgsCurve * interiorRing(int i) const
Retrieves an interior ring from the curve polygon.
Abstract base class for curved geometry type.
Definition qgscurve.h:36
virtual int numPoints() const =0
Returns the number of points in the curve.
QgsCurve * segmentize(double tolerance=M_PI_2/90, SegmentationToleranceType toleranceType=MaximumAngle) const override
Returns a geometry without curves.
Definition qgscurve.cpp:175
virtual bool isClosed() const
Returns true if the curve is closed.
Definition qgscurve.cpp:53
QgsCurve * clone() const override=0
Clones the geometry by performing a deep copy.
virtual QgsPolygon * toPolygon(unsigned int segments=36) const
Returns a segmented polygon.
Base class for feedback objects to be used for cancellation of something running in a worker thread.
Definition qgsfeedback.h:45
bool isCanceled() const
Tells whether the operation has been canceled already.
Definition qgsfeedback.h:54
int partCount() const override
Returns count of parts contained in the geometry.
int numGeometries() const
Returns the number of geometries within the collection.
const QgsAbstractGeometry * geometryN(int n) const
Returns a const reference to a geometry from within the collection.
static QgsPoint pointOnLineWithDistance(const QgsPoint &startPoint, const QgsPoint &directionPoint, double distance)
Returns a point a specified distance toward a second point.
static double triangleArea(double aX, double aY, double bX, double bY, double cX, double cY)
Returns the area of the triangle denoted by the points (aX, aY), (bX, bY) and (cX,...
static void weightedPointInTriangle(double aX, double aY, double bX, double bY, double cX, double cY, double weightB, double weightC, double &pointX, double &pointY)
Returns a weighted point inside the triangle denoted by the points (aX, aY), (bX, bY) and (cX,...
static double lineAngle(double x1, double y1, double x2, double y2)
Calculates the direction of line joining two points in radians, clockwise from the north direction.
static bool pointContinuesArc(const QgsPoint &a1, const QgsPoint &a2, const QgsPoint &a3, const QgsPoint &b, double distanceTolerance, double pointSpacingAngleTolerance)
Returns true if point b is on the arc formed by points a1, a2, and a3, but not within that arc portio...
static int circleCircleOuterTangents(const QgsPointXY &center1, double radius1, const QgsPointXY &center2, double radius2, QgsPointXY &line1P1, QgsPointXY &line1P2, QgsPointXY &line2P1, QgsPointXY &line2P2)
Calculates the outer tangent points for two circles, centered at center1 and center2 and with radii o...
static int leftOfLine(const double x, const double y, const double x1, const double y1, const double x2, const double y2)
Returns a value < 0 if the point (x, y) is left of the line from (x1, y1) -> (x2, y2).
static void circleCenterRadius(const QgsPoint &pt1, const QgsPoint &pt2, const QgsPoint &pt3, double &radius, double &centerX, double &centerY)
Returns radius and center of the circle through pt1, pt2, pt3.
A geometry is the spatial representation of a feature.
static QgsGeometry fromRect(const QgsRectangle &rect)
Creates a new geometry from a QgsRectangle.
static QgsGeometry collectGeometry(const QVector< QgsGeometry > &geometries)
Creates a new multipart geometry from a list of QgsGeometry objects.
QgsPointXY closestVertex(const QgsPointXY &point, int &closestVertexIndex, int &previousVertexIndex, int &nextVertexIndex, double &sqrDist) const
Returns the vertex closest to the given point, the corresponding vertex index, squared distance snap ...
const QgsAbstractGeometry * constGet() const
Returns a non-modifiable (const) reference to the underlying abstract geometry primitive.
Qgis::GeometryType type
bool isMultipart() const
Returns true if WKB of the geometry is of WKBMulti* type.
static QgsGeometryEngine * createGeometryEngine(const QgsAbstractGeometry *geometry)
Creates and returns a new geometry engine representing the specified geometry.
Qgis::GeometryOperationResult addPart(const QVector< QgsPointXY > &points, Qgis::GeometryType geomType=Qgis::GeometryType::Unknown)
Adds a new part to a the geometry.
bool removeDuplicateNodes(double epsilon=4 *std::numeric_limits< double >::epsilon(), bool useZValues=false)
Removes duplicate nodes from the geometry, wherever removing the nodes does not result in a degenerat...
static QgsGeometry unaryUnion(const QVector< QgsGeometry > &geometries, const QgsGeometryParameters &parameters=QgsGeometryParameters())
Compute the unary union on a list of geometries.
QgsRectangle boundingBox() const
Returns the bounding box of the geometry.
Qgis::GeometryOperationResult rotate(double rotation, const QgsPointXY &center)
Rotate this geometry around the Z axis.
QgsGeometry triangularWavesRandomized(double minimumWavelength, double maximumWavelength, double minimumAmplitude, double maximumAmplitude, unsigned long seed=0) const
Constructs randomized triangular waves along the boundary of the geometry, with the specified wavelen...
QgsGeometry triangularWaves(double wavelength, double amplitude, bool strictWavelength=false) const
Constructs triangular waves along the boundary of the geometry, with the specified wavelength and amp...
QgsInternalGeometryEngine(const QgsGeometry &geometry)
The caller is responsible that the geometry is available and unchanged for the whole lifetime of this...
QgsGeometry roundWaves(double wavelength, double amplitude, bool strictWavelength=false) const
Constructs rounded (sine-like) waves along the boundary of the geometry, with the specified wavelengt...
QgsGeometry poleOfInaccessibility(double precision, double *distanceFromBoundary=nullptr) const
Calculates the approximate pole of inaccessibility for a surface, which is the most distant internal ...
QgsGeometry squareWaves(double wavelength, double amplitude, bool strictWavelength=false) const
Constructs square waves along the boundary of the geometry, with the specified wavelength and amplitu...
QgsGeometry variableWidthBufferByM(int segments) const
Calculates a variable width buffer using the m-values from a (multi)line geometry.
QgsGeometry extrude(double x, double y) const
Will extrude a line or (segmentized) curve by a given offset and return a polygon representation of i...
QgsGeometry roundWavesRandomized(double minimumWavelength, double maximumWavelength, double minimumAmplitude, double maximumAmplitude, unsigned long seed=0) const
Constructs randomized rounded (sine-like) waves along the boundary of the geometry,...
QgsGeometry orthogonalize(double tolerance=1.0E-8, int maxIterations=1000, double angleThreshold=15.0) const
Attempts to orthogonalize a line or polygon geometry by shifting vertices to make the geometries angl...
QgsGeometry variableWidthBuffer(int segments, const std::function< std::unique_ptr< double[] >(const QgsLineString *line) > &widthFunction) const
Calculates a variable width buffer for a (multi)curve geometry.
QString lastError() const
Returns an error string referring to the last error encountered.
QgsGeometry orientedMinimumBoundingBox(double &area, double &angle, double &width, double &height) const
Returns the oriented minimum bounding box for the geometry, which is the smallest (by area) rotated r...
QgsGeometry densifyByDistance(double distance) const
Densifies the geometry by adding regularly placed extra nodes inside each segment so that the maximum...
QgsGeometry taperedBuffer(double startWidth, double endWidth, int segments) const
Calculates a tapered width buffer for a (multi)curve geometry.
QgsGeometry densifyByCount(int extraNodesPerSegment) const
Densifies the geometry by adding the specified number of extra nodes within each segment of the geome...
QgsGeometry applyDashPattern(const QVector< double > &pattern, Qgis::DashPatternLineEndingRule startRule=Qgis::DashPatternLineEndingRule::NoRule, Qgis::DashPatternLineEndingRule endRule=Qgis::DashPatternLineEndingRule::NoRule, Qgis::DashPatternSizeAdjustment adjustment=Qgis::DashPatternSizeAdjustment::ScaleBothDashAndGap, double patternOffset=0) const
Applies a dash pattern to a geometry, returning a MultiLineString geometry which is the input geometr...
static QVector< QgsPointXY > randomPointsInPolygon(const QgsGeometry &polygon, int count, const std::function< bool(const QgsPointXY &) > &acceptPoint, unsigned long seed=0, QgsFeedback *feedback=nullptr, int maxTriesPerPoint=0)
Returns a list of count random points generated inside a polygon geometry (if acceptPoint is specifie...
QgsGeometry squareWavesRandomized(double minimumWavelength, double maximumWavelength, double minimumAmplitude, double maximumAmplitude, unsigned long seed=0) const
Constructs randomized square waves along the boundary of the geometry, with the specified wavelength ...
QgsGeometry convertToCurves(double distanceTolerance, double angleTolerance) const
Attempts to convert a non-curved geometry into a curved geometry type (e.g.
bool isAxisParallelRectangle(double maximumDeviation, bool simpleRectanglesOnly=false) const
Returns true if the geometry is a polygon that is almost an axis-parallel rectangle.
Represents a single 2D line segment, consisting of a 2D start and end vertex only.
int pointLeftOfLine(const QgsPointXY &point) const
Tests if a point is to the left of the line segment.
double endY() const
Returns the segment's end y-coordinate.
QgsPointXY end() const
Returns the segment's end point.
double endX() const
Returns the segment's end x-coordinate.
QgsPointXY start() const
Returns the segment's start point.
double startX() const
Returns the segment's start x-coordinate.
void reverse()
Reverses the line segment, so that the start and end points are flipped.
double startY() const
Returns the segment's start y-coordinate.
Line string geometry type, with support for z-dimension and m-values.
bool isClosed() const override
Returns true if the curve is closed.
const double * yData() const
Returns a const pointer to the y vertex data.
const double * xData() const
Returns a const pointer to the x vertex data.
double length() const override
Returns the planar, 2-dimensional length of the geometry.
QgsPoint startPoint() const override
Returns the starting point of the curve.
static QgsLineString * fromBezierCurve(const QgsPoint &start, const QgsPoint &controlPoint1, const QgsPoint &controlPoint2, const QgsPoint &end, int segments=30)
Returns a new linestring created by segmentizing the bezier curve between start and end,...
int numPoints() const override
Returns the number of points in the curve.
QgsPoint pointN(int i) const
Returns the specified point from inside the line string.
int nCoordinates() const override
Returns the number of nodes contained in the geometry.
double yAt(int index) const override
Returns the y-coordinate of the specified node in the line string.
void setYAt(int index, double y)
Sets the y-coordinate of the specified node in the line string.
double mAt(int index) const override
Returns the m value of the specified node in the line string.
void setXAt(int index, double x)
Sets the x-coordinate of the specified node in the line string.
QgsLineString * clone() const override
Clones the geometry by performing a deep copy.
double zAt(int index) const override
Returns the z-coordinate of the specified node in the line string.
double xAt(int index) const override
Returns the x-coordinate of the specified node in the line string.
Multi curve geometry collection.
Multi polygon geometry collection.
Multi surface geometry collection.
A class to represent a 2D point.
Definition qgspointxy.h:59
double y
Definition qgspointxy.h:63
double x
Definition qgspointxy.h:62
Point geometry type, with support for z-dimension and m-values.
Definition qgspoint.h:49
QgsPoint * clone() const override
Clones the geometry by performing a deep copy.
Definition qgspoint.cpp:105
QgsPoint vertexAt(QgsVertexId) const override
Returns the point corresponding to a specified vertex id.
Definition qgspoint.cpp:526
double x
Definition qgspoint.h:52
bool isEmpty() const override
Returns true if the geometry is empty.
Definition qgspoint.cpp:773
double distance(double x, double y) const
Returns the Cartesian 2D distance between this point and a specified x, y coordinate.
Definition qgspoint.h:343
double y
Definition qgspoint.h:53
bool nextVertex(QgsVertexId &id, QgsPoint &vertex) const override
Returns next vertex id and coordinates.
Definition qgspoint.cpp:477
Polygon geometry type.
Definition qgspolygon.h:34
void setExteriorRing(QgsCurve *ring) override
Sets the exterior ring of the polygon.
void addInteriorRing(QgsCurve *ring) override
Adds an interior ring to the geometry (takes ownership)
bool intersects(const QgsLineSegment2D &segment, QgsPointXY &intersectPoint) const
Finds the closest intersection point of the ray and a line segment.
A rectangle specified with double values.
double xMinimum() const
Returns the x minimum value (left side of rectangle).
double yMinimum() const
Returns the y minimum value (bottom side of rectangle).
double width() const
Returns the width of the rectangle.
double xMaximum() const
Returns the x maximum value (right side of rectangle).
double yMaximum() const
Returns the y maximum value (top side of rectangle).
double height() const
Returns the height of the rectangle.
Surface geometry type.
Definition qgssurface.h:34
Class that takes care of tessellation of polygons into triangles.
QVector< float > data() const
Returns array of triangle vertex data.
void addPolygon(const QgsPolygon &polygon, float extrusionHeight)
Tessellates a triangle and adds its vertex entries to the output data array.
int dataVerticesCount() const
Returns the number of vertices stored in the output data array.
A class to represent a vector.
Definition qgsvector.h:30
double lengthSquared() const
Returns the length of the vector.
Definition qgsvector.h:138
double crossProduct(QgsVector v) const
Returns the 2D cross product of this vector and another vector v.
Definition qgsvector.h:192
QgsVector normalized() const
Returns the vector's normalized (or "unit") vector (ie same angle but length of 1....
Definition qgsvector.cpp:28
double length() const
Returns the length of the vector.
Definition qgsvector.h:128
static Qgis::GeometryType geometryType(Qgis::WkbType type)
Returns the geometry type for a WKB type, e.g., both MultiPolygon and CurvePolygon would have a Polyg...
static bool isCurvedType(Qgis::WkbType type)
Returns true if the WKB type is a curved type or can contain curved geometries.
static Qgis::WkbType flatType(Qgis::WkbType type)
Returns the flat type for a WKB type.
As part of the API refactoring and improvements which landed in the Processing API was substantially reworked from the x version This was done in order to allow much of the underlying Processing framework to be ported into c
bool qgsDoubleNear(double a, double b, double epsilon=4 *std::numeric_limits< double >::epsilon())
Compare two doubles (but allow some difference)
Definition qgis.h:4332
std::unique_ptr< QgsCurve > lineToCurve(const QgsCurve *curve, double distanceTolerance, double pointSpacingAngleTolerance)
bool matchesOrientation(std::array< Direction, 4 > dirs, std::array< Direction, 4 > oriented)
double normalizedDotProduct(const QgsPoint &a, const QgsPoint &b, const QgsPoint &c)
bool isClockwise(std::array< Direction, 4 > dirs)
Checks whether the 4 directions in dirs make up a clockwise rectangle.
std::unique_ptr< QgsLineString > squareWavesRandomizedAlongLine(const QgsLineString *line, const double minimumWavelength, const double maximumWavelength, const double minimumAmplitude, const double maximumAmplitude, std::uniform_real_distribution<> &uniformDist, std::mt19937 &mt)
std::unique_ptr< QgsLineString > squareWavesAlongLine(const QgsLineString *line, double wavelength, const double amplitude, const bool strictWavelength)
QgsAbstractGeometry * orthogonalizeGeom(const QgsAbstractGeometry *geom, int maxIterations, double tolerance, double lowerThreshold, double upperThreshold)
std::unique_ptr< QgsLineString > triangularWavesRandomizedAlongLine(const QgsLineString *line, const double minimumWavelength, const double maximumWavelength, const double minimumAmplitude, const double maximumAmplitude, std::uniform_real_distribution<> &uniformDist, std::mt19937 &mt)
QgsLineString * doDensify(const QgsLineString *ring, int extraNodesPerSegment=-1, double distance=1)
Direction getEdgeDirection(const QgsPoint &p1, const QgsPoint &p2, double maxDev)
Determines the direction of an edge from p1 to p2.
std::unique_ptr< QgsAbstractGeometry > applyDashPatternPrivate(const QgsAbstractGeometry *geom, const QVector< double > &pattern, Qgis::DashPatternLineEndingRule startRule, Qgis::DashPatternLineEndingRule endRule, Qgis::DashPatternSizeAdjustment adjustment, double patternOffset)
std::unique_ptr< QgsLineString > roundWavesRandomizedAlongLine(const QgsLineString *line, const double minimumWavelength, const double maximumWavelength, const double minimumAmplitude, const double maximumAmplitude, std::uniform_real_distribution<> &uniformDist, std::mt19937 &mt)
QgsVector calcMotion(const QgsPoint &a, const QgsPoint &b, const QgsPoint &c, double lowerThreshold, double upperThreshold)
QVector< QgsPointXY > generateSegmentCurve(const QgsPoint &center1, const double radius1, const QgsPoint &center2, const double radius2)
bool dotProductWithinAngleTolerance(double dotProduct, double lowerThreshold, double upperThreshold)
QgsLineString * doOrthogonalize(QgsLineString *ring, int iterations, double tolerance, double lowerThreshold, double upperThreshold)
double squareness(QgsLineString *ring, double lowerThreshold, double upperThreshold)
std::unique_ptr< QgsMultiLineString > dashPatternAlongLine(const QgsLineString *line, const QVector< double > &pattern, Qgis::DashPatternLineEndingRule startRule, Qgis::DashPatternLineEndingRule endRule, Qgis::DashPatternSizeAdjustment adjustment, double patternOffset)
std::unique_ptr< QgsAbstractGeometry > convertGeometryToCurves(const QgsAbstractGeometry *geom, double distanceTolerance, double angleTolerance)
std::unique_ptr< QgsAbstractGeometry > roundWavesPrivate(const QgsAbstractGeometry *geom, double wavelength, double amplitude, bool strictWavelength)
bool isCounterClockwise(std::array< Direction, 4 > dirs)
Checks whether the 4 directions in dirs make up a counter-clockwise rectangle.
QgsAbstractGeometry * densifyGeometry(const QgsAbstractGeometry *geom, int extraNodesPerSegment=1, double distance=1)
std::unique_ptr< QgsLineString > roundWavesAlongLine(const QgsLineString *line, double wavelength, const double amplitude, const bool strictWavelength)
std::pair< bool, std::array< Direction, 4 > > getEdgeDirections(const QgsPolygon *g, double maxDev)
Checks whether the polygon consists of four nearly axis-parallel sides.
std::unique_ptr< QgsAbstractGeometry > squareWavesRandomizedPrivate(const QgsAbstractGeometry *geom, double minimumWavelength, double maximumWavelength, double minimumAmplitude, double maximumAmplitude, std::uniform_real_distribution<> &uniformDist, std::mt19937 &mt)
std::unique_ptr< QgsAbstractGeometry > triangularWavesRandomizedPrivate(const QgsAbstractGeometry *geom, double minimumWavelength, double maximumWavelength, double minimumAmplitude, double maximumAmplitude, std::uniform_real_distribution<> &uniformDist, std::mt19937 &mt)
std::unique_ptr< QgsAbstractGeometry > triangularWavesPrivate(const QgsAbstractGeometry *geom, double wavelength, double amplitude, bool strictWavelength)
std::unique_ptr< QgsAbstractGeometry > squareWavesPrivate(const QgsAbstractGeometry *geom, double wavelength, double amplitude, bool strictWavelength)
std::unique_ptr< QgsLineString > triangularWavesAlongLine(const QgsLineString *line, double wavelength, const double amplitude, const bool strictWavelength)
std::unique_ptr< QgsAbstractGeometry > roundWavesRandomizedPrivate(const QgsAbstractGeometry *geom, double minimumWavelength, double maximumWavelength, double minimumAmplitude, double maximumAmplitude, std::uniform_real_distribution<> &uniformDist, std::mt19937 &mt)
QLineF segment(int index, QRectF rect, double radius)
int precision
Utility class for identifying a unique vertex within a geometry.
Definition qgsvertexid.h:31