QGIS API Documentation  2.18.21-Las Palmas (9fba24a)
qgscircularstringv2.cpp
Go to the documentation of this file.
1 /***************************************************************************
2  qgscircularstringv2.cpp
3  -----------------------
4  begin : September 2014
5  copyright : (C) 2014 by Marco Hugentobler
6  email : marco at sourcepole dot ch
7  ***************************************************************************/
8 
9 /***************************************************************************
10  * *
11  * This program is free software; you can redistribute it and/or modify *
12  * it under the terms of the GNU General Public License as published by *
13  * the Free Software Foundation; either version 2 of the License, or *
14  * (at your option) any later version. *
15  * *
16  ***************************************************************************/
17 
18 #include "qgscircularstringv2.h"
19 #include "qgsapplication.h"
20 #include "qgscoordinatetransform.h"
21 #include "qgsgeometryutils.h"
22 #include "qgslinestringv2.h"
23 #include "qgsmaptopixel.h"
24 #include "qgspointv2.h"
25 #include "qgswkbptr.h"
26 #include <QPainter>
27 #include <QPainterPath>
28 
30 {
32 }
33 
35 {
36 
37 }
38 
39 bool QgsCircularStringV2::operator==( const QgsCurveV2& other ) const
40 {
41  const QgsCircularStringV2* otherLine = dynamic_cast< const QgsCircularStringV2* >( &other );
42  if ( !otherLine )
43  return false;
44 
45  return *otherLine == *this;
46 }
47 
48 bool QgsCircularStringV2::operator!=( const QgsCurveV2& other ) const
49 {
50  return !operator==( other );
51 }
52 
54 {
55  return new QgsCircularStringV2( *this );
56 }
57 
59 {
61  mX.clear();
62  mY.clear();
63  mZ.clear();
64  mM.clear();
65  clearCache();
66 }
67 
69 {
70  QgsRectangle bbox;
71  int nPoints = numPoints();
72  for ( int i = 0; i < ( nPoints - 2 ) ; i += 2 )
73  {
74  if ( i == 0 )
75  {
76  bbox = segmentBoundingBox( QgsPointV2( mX[i], mY[i] ), QgsPointV2( mX[i + 1], mY[i + 1] ), QgsPointV2( mX[i + 2], mY[i + 2] ) );
77  }
78  else
79  {
80  QgsRectangle segmentBox = segmentBoundingBox( QgsPointV2( mX[i], mY[i] ), QgsPointV2( mX[i + 1], mY[i + 1] ), QgsPointV2( mX[i + 2], mY[i + 2] ) );
81  bbox.combineExtentWith( segmentBox );
82  }
83  }
84 
85  if ( nPoints > 0 && nPoints % 2 == 0 )
86  {
87  if ( nPoints == 2 )
88  {
89  bbox.combineExtentWith( mX[ 0 ], mY[ 0 ] );
90  }
91  bbox.combineExtentWith( mX[ nPoints - 1 ], mY[ nPoints - 1 ] );
92  }
93  return bbox;
94 }
95 
96 QgsRectangle QgsCircularStringV2::segmentBoundingBox( const QgsPointV2& pt1, const QgsPointV2& pt2, const QgsPointV2& pt3 )
97 {
98  double centerX, centerY, radius;
99  QgsGeometryUtils::circleCenterRadius( pt1, pt2, pt3, radius, centerX, centerY );
100 
101  double p1Angle = QgsGeometryUtils::ccwAngle( pt1.y() - centerY, pt1.x() - centerX );
102  if ( p1Angle > 360 )
103  {
104  p1Angle -= 360;
105  }
106  double p2Angle = QgsGeometryUtils::ccwAngle( pt2.y() - centerY, pt2.x() - centerX );
107  if ( p2Angle > 360 )
108  {
109  p2Angle -= 360;
110  }
111  double p3Angle = QgsGeometryUtils::ccwAngle( pt3.y() - centerY, pt3.x() - centerX );
112  if ( p3Angle > 360 )
113  {
114  p3Angle -= 360;
115  }
116 
117  //start point, end point and compass points in between can be on bounding box
118  QgsRectangle bbox( pt1.x(), pt1.y(), pt1.x(), pt1.y() );
119  bbox.combineExtentWith( pt3.x(), pt3.y() );
120 
121  QgsPointSequenceV2 compassPoints = compassPointsOnSegment( p1Angle, p2Angle, p3Angle, centerX, centerY, radius );
122  QgsPointSequenceV2::const_iterator cpIt = compassPoints.constBegin();
123  for ( ; cpIt != compassPoints.constEnd(); ++cpIt )
124  {
125  bbox.combineExtentWith( cpIt->x(), cpIt->y() );
126  }
127  return bbox;
128 }
129 
130 QgsPointSequenceV2 QgsCircularStringV2::compassPointsOnSegment( double p1Angle, double p2Angle, double p3Angle, double centerX, double centerY, double radius )
131 {
132  QgsPointSequenceV2 pointList;
133 
134  QgsPointV2 nPoint( centerX, centerY + radius );
135  QgsPointV2 ePoint( centerX + radius, centerY );
136  QgsPointV2 sPoint( centerX, centerY - radius );
137  QgsPointV2 wPoint( centerX - radius, centerY );
138 
139  if ( p3Angle >= p1Angle )
140  {
141  if ( p2Angle > p1Angle && p2Angle < p3Angle )
142  {
143  if ( p1Angle <= 90 && p3Angle >= 90 )
144  {
145  pointList.append( nPoint );
146  }
147  if ( p1Angle <= 180 && p3Angle >= 180 )
148  {
149  pointList.append( wPoint );
150  }
151  if ( p1Angle <= 270 && p3Angle >= 270 )
152  {
153  pointList.append( sPoint );
154  }
155  }
156  else
157  {
158  pointList.append( ePoint );
159  if ( p1Angle >= 90 || p3Angle <= 90 )
160  {
161  pointList.append( nPoint );
162  }
163  if ( p1Angle >= 180 || p3Angle <= 180 )
164  {
165  pointList.append( wPoint );
166  }
167  if ( p1Angle >= 270 || p3Angle <= 270 )
168  {
169  pointList.append( sPoint );
170  }
171  }
172  }
173  else
174  {
175  if ( p2Angle < p1Angle && p2Angle > p3Angle )
176  {
177  if ( p1Angle >= 270 && p3Angle <= 270 )
178  {
179  pointList.append( sPoint );
180  }
181  if ( p1Angle >= 180 && p3Angle <= 180 )
182  {
183  pointList.append( wPoint );
184  }
185  if ( p1Angle >= 90 && p3Angle <= 90 )
186  {
187  pointList.append( nPoint );
188  }
189  }
190  else
191  {
192  pointList.append( ePoint );
193  if ( p1Angle <= 270 || p3Angle >= 270 )
194  {
195  pointList.append( sPoint );
196  }
197  if ( p1Angle <= 180 || p3Angle >= 180 )
198  {
199  pointList.append( wPoint );
200  }
201  if ( p1Angle <= 90 || p3Angle >= 90 )
202  {
203  pointList.append( nPoint );
204  }
205  }
206  }
207  return pointList;
208 }
209 
211 {
212  if ( !wkbPtr )
213  return false;
214 
215  QgsWKBTypes::Type type = wkbPtr.readHeader();
217  {
218  return false;
219  }
220  clearCache();
221  mWkbType = type;
222 
223  //type
224  bool hasZ = is3D();
225  bool hasM = isMeasure();
226  int nVertices = 0;
227  wkbPtr >> nVertices;
228  mX.resize( nVertices );
229  mY.resize( nVertices );
230  hasZ ? mZ.resize( nVertices ) : mZ.clear();
231  hasM ? mM.resize( nVertices ) : mM.clear();
232  for ( int i = 0; i < nVertices; ++i )
233  {
234  wkbPtr >> mX[i];
235  wkbPtr >> mY[i];
236  if ( hasZ )
237  {
238  wkbPtr >> mZ[i];
239  }
240  if ( hasM )
241  {
242  wkbPtr >> mM[i];
243  }
244  }
245 
246  return true;
247 }
248 
250 {
251  clear();
252 
254 
255  if ( QgsWKBTypes::flatType( parts.first ) != QgsWKBTypes::CircularString )
256  return false;
257  mWkbType = parts.first;
258 
259  setPoints( QgsGeometryUtils::pointsFromWKT( parts.second, is3D(), isMeasure() ) );
260  return true;
261 }
262 
264 {
265  int size = sizeof( char ) + sizeof( quint32 ) + sizeof( quint32 );
266  size += numPoints() * ( 2 + is3D() + isMeasure() ) * sizeof( double );
267  return size;
268 }
269 
270 unsigned char* QgsCircularStringV2::asWkb( int& binarySize ) const
271 {
272  binarySize = wkbSize();
273  unsigned char* geomPtr = new unsigned char[binarySize];
274  QgsWkbPtr wkb( geomPtr, binarySize );
275  wkb << static_cast<char>( QgsApplication::endian() );
276  wkb << static_cast<quint32>( wkbType() );
277  QgsPointSequenceV2 pts;
278  points( pts );
279  QgsGeometryUtils::pointsToWKB( wkb, pts, is3D(), isMeasure() );
280  return geomPtr;
281 }
282 
283 QString QgsCircularStringV2::asWkt( int precision ) const
284 {
285  QString wkt = wktTypeStr() + ' ';
286  QgsPointSequenceV2 pts;
287  points( pts );
288  wkt += QgsGeometryUtils::pointsToWKT( pts, precision, is3D(), isMeasure() );
289  return wkt;
290 }
291 
292 QDomElement QgsCircularStringV2::asGML2( QDomDocument& doc, int precision, const QString& ns ) const
293 {
294  // GML2 does not support curves
295  QgsLineStringV2* line = curveToLine();
296  QDomElement gml = line->asGML2( doc, precision, ns );
297  delete line;
298  return gml;
299 }
300 
301 QDomElement QgsCircularStringV2::asGML3( QDomDocument& doc, int precision, const QString& ns ) const
302 {
303  QgsPointSequenceV2 pts;
304  points( pts );
305 
306  QDomElement elemCurve = doc.createElementNS( ns, "Curve" );
307  QDomElement elemSegments = doc.createElementNS( ns, "segments" );
308  QDomElement elemArcString = doc.createElementNS( ns, "ArcString" );
309  elemArcString.appendChild( QgsGeometryUtils::pointsToGML3( pts, doc, precision, ns, is3D() ) );
310  elemSegments.appendChild( elemArcString );
311  elemCurve.appendChild( elemSegments );
312  return elemCurve;
313 }
314 
315 QString QgsCircularStringV2::asJSON( int precision ) const
316 {
317  // GeoJSON does not support curves
318  QgsLineStringV2* line = curveToLine();
319  QString json = line->asJSON( precision );
320  delete line;
321  return json;
322 }
323 
324 //curve interface
326 {
327  int nPoints = numPoints();
328  double length = 0;
329  for ( int i = 0; i < ( nPoints - 2 ) ; i += 2 )
330  {
331  length += QgsGeometryUtils::circleLength( mX[i], mY[i], mX[i + 1], mY[i + 1], mX[i + 2], mY[i + 2] );
332  }
333  return length;
334 }
335 
337 {
338  if ( numPoints() < 1 )
339  {
340  return QgsPointV2();
341  }
342  return pointN( 0 );
343 }
344 
346 {
347  if ( numPoints() < 1 )
348  {
349  return QgsPointV2();
350  }
351  return pointN( numPoints() - 1 );
352 }
353 
355 {
356  QgsLineStringV2* line = new QgsLineStringV2();
358  int nPoints = numPoints();
359 
360  for ( int i = 0; i < ( nPoints - 2 ) ; i += 2 )
361  {
362  segmentize( pointN( i ), pointN( i + 1 ), pointN( i + 2 ), points, tolerance, toleranceType );
363  }
364 
365  line->setPoints( points );
366  return line;
367 }
368 
370 {
371  return qMin( mX.size(), mY.size() );
372 }
373 
375 {
376  if ( qMin( mX.size(), mY.size() ) <= i )
377  {
378  return QgsPointV2();
379  }
380 
381  double x = mX.at( i );
382  double y = mY.at( i );
383  double z = 0;
384  double m = 0;
385 
386  if ( is3D() )
387  {
388  z = mZ.at( i );
389  }
390  if ( isMeasure() )
391  {
392  m = mM.at( i );
393  }
394 
396  if ( is3D() && isMeasure() )
397  {
399  }
400  else if ( is3D() )
401  {
403  }
404  else if ( isMeasure() )
405  {
407  }
408  return QgsPointV2( t, x, y, z, m );
409 }
410 
412 {
413  pts.clear();
414  int nPts = numPoints();
415  for ( int i = 0; i < nPts; ++i )
416  {
417  pts.push_back( pointN( i ) );
418  }
419 }
420 
422 {
423  clearCache();
424 
425  if ( points.size() < 1 )
426  {
428  mX.clear();
429  mY.clear();
430  mZ.clear();
431  mM.clear();
432  return;
433  }
434 
435  //get wkb type from first point
436  const QgsPointV2& firstPt = points.at( 0 );
437  bool hasZ = firstPt.is3D();
438  bool hasM = firstPt.isMeasure();
439 
441 
442  mX.resize( points.size() );
443  mY.resize( points.size() );
444  if ( hasZ )
445  {
446  mZ.resize( points.size() );
447  }
448  else
449  {
450  mZ.clear();
451  }
452  if ( hasM )
453  {
454  mM.resize( points.size() );
455  }
456  else
457  {
458  mM.clear();
459  }
460 
461  for ( int i = 0; i < points.size(); ++i )
462  {
463  mX[i] = points[i].x();
464  mY[i] = points[i].y();
465  if ( hasZ )
466  {
467  mZ[i] = points[i].z();
468  }
469  if ( hasM )
470  {
471  mM[i] = points[i].m();
472  }
473  }
474 }
475 
476 void QgsCircularStringV2::segmentize( const QgsPointV2& p1, const QgsPointV2& p2, const QgsPointV2& p3, QgsPointSequenceV2 &points, double tolerance, SegmentationToleranceType toleranceType ) const
477 {
478  bool clockwise = false;
479  int segSide = segmentSide( p1, p3, p2 );
480  if ( segSide == -1 )
481  {
482  clockwise = true;
483  }
484 
485  QgsPointV2 circlePoint1 = clockwise ? p3 : p1;
486  QgsPointV2 circlePoint2 = p2;
487  QgsPointV2 circlePoint3 = clockwise ? p1 : p3 ;
488 
489  //adapted code from postgis
490  double radius = 0;
491  double centerX = 0;
492  double centerY = 0;
493  QgsGeometryUtils::circleCenterRadius( circlePoint1, circlePoint2, circlePoint3, radius, centerX, centerY );
494 
495 
496  if ( circlePoint1 != circlePoint3 && ( radius < 0 || qgsDoubleNear( segSide, 0.0 ) ) ) //points are colinear
497  {
498  points.append( p1 );
499  points.append( p2 );
500  points.append( p3 );
501  return;
502  }
503 
504  double increment = tolerance; //one segment per degree
505  if ( toleranceType == QgsAbstractGeometryV2::MaximumDifference )
506  {
507  double halfAngle = acos( -tolerance / radius + 1 );
508  increment = 2 * halfAngle;
509  }
510 
511  //angles of pt1, pt2, pt3
512  double a1 = atan2( circlePoint1.y() - centerY, circlePoint1.x() - centerX );
513  double a2 = atan2( circlePoint2.y() - centerY, circlePoint2.x() - centerX );
514  double a3 = atan2( circlePoint3.y() - centerY, circlePoint3.x() - centerX );
515 
516  /* Adjust a3 up so we can increment from a1 to a3 cleanly */
517  if ( a3 <= a1 )
518  a3 += 2.0 * M_PI;
519  if ( a2 < a1 )
520  a2 += 2.0 * M_PI;
521 
522  bool hasZ = is3D();
523  bool hasM = isMeasure();
524 
525  double x, y;
526  double z = 0;
527  double m = 0;
528 
529  QList<QgsPointV2> stringPoints;
530  stringPoints.insert( clockwise ? 0 : stringPoints.size(), circlePoint1 );
531  if ( circlePoint2 != circlePoint3 && circlePoint1 != circlePoint2 ) //draw straight line segment if two points have the same position
532  {
533  QgsWKBTypes::Type pointWkbType = QgsWKBTypes::Point;
534  if ( hasZ )
535  pointWkbType = QgsWKBTypes::addZ( pointWkbType );
536  if ( hasM )
537  pointWkbType = QgsWKBTypes::addM( pointWkbType );
538 
539  //make sure the curve point p2 is part of the segmentized vertices. But only if p1 != p3
540  bool addP2 = true;
541  if ( qgsDoubleNear( circlePoint1.x(), circlePoint3.x() ) && qgsDoubleNear( circlePoint1.y(), circlePoint3.y() ) )
542  {
543  addP2 = false;
544  }
545 
546  for ( double angle = a1 + increment; angle < a3; angle += increment )
547  {
548  if (( addP2 && angle > a2 ) )
549  {
550  stringPoints.insert( clockwise ? 0 : stringPoints.size(), circlePoint2 );
551  addP2 = false;
552  }
553 
554  x = centerX + radius * cos( angle );
555  y = centerY + radius * sin( angle );
556 
557  if ( !hasZ && !hasM )
558  {
559  stringPoints.insert( clockwise ? 0 : stringPoints.size(), QgsPointV2( x, y ) );
560  continue;
561  }
562 
563  if ( hasZ )
564  {
565  z = interpolateArc( angle, a1, a2, a3, circlePoint1.z(), circlePoint2.z(), circlePoint3.z() );
566  }
567  if ( hasM )
568  {
569  m = interpolateArc( angle, a1, a2, a3, circlePoint1.m(), circlePoint2.m(), circlePoint3.m() );
570  }
571 
572  stringPoints.insert( clockwise ? 0 : stringPoints.size(), QgsPointV2( pointWkbType, x, y, z, m ) );
573  }
574  }
575  stringPoints.insert( clockwise ? 0 : stringPoints.size(), circlePoint3 );
576  points.append( stringPoints );
577 }
578 
579 int QgsCircularStringV2::segmentSide( const QgsPointV2& pt1, const QgsPointV2& pt3, const QgsPointV2& pt2 ) const
580 {
581  double side = (( pt2.x() - pt1.x() ) * ( pt3.y() - pt1.y() ) - ( pt3.x() - pt1.x() ) * ( pt2.y() - pt1.y() ) );
582  if ( side == 0.0 )
583  {
584  return 0;
585  }
586  else
587  {
588  if ( side < 0 )
589  {
590  return -1;
591  }
592  if ( side > 0 )
593  {
594  return 1;
595  }
596  return 0;
597  }
598 }
599 
600 double QgsCircularStringV2::interpolateArc( double angle, double a1, double a2, double a3, double zm1, double zm2, double zm3 ) const
601 {
602  /* Counter-clockwise sweep */
603  if ( a1 < a2 )
604  {
605  if ( angle <= a2 )
606  return zm1 + ( zm2 - zm1 ) * ( angle - a1 ) / ( a2 - a1 );
607  else
608  return zm2 + ( zm3 - zm2 ) * ( angle - a2 ) / ( a3 - a2 );
609  }
610  /* Clockwise sweep */
611  else
612  {
613  if ( angle >= a2 )
614  return zm1 + ( zm2 - zm1 ) * ( a1 - angle ) / ( a1 - a2 );
615  else
616  return zm2 + ( zm3 - zm2 ) * ( a2 - angle ) / ( a2 - a3 );
617  }
618 }
619 
621 {
622  QPainterPath path;
623  addToPainterPath( path );
624  p.drawPath( path );
625 }
626 
628 {
629  clearCache();
630 
631  double* zArray = mZ.data();
632 
633  bool hasZ = is3D();
634  int nPoints = numPoints();
635  bool useDummyZ = !hasZ || !transformZ;
636  if ( useDummyZ )
637  {
638  zArray = new double[nPoints];
639  for ( int i = 0; i < nPoints; ++i )
640  {
641  zArray[i] = 0;
642  }
643  }
644  ct.transformCoords( nPoints, mX.data(), mY.data(), zArray, d );
645  if ( useDummyZ )
646  {
647  delete[] zArray;
648  }
649 }
650 
652 {
653  clearCache();
654 
655  int nPoints = numPoints();
656  for ( int i = 0; i < nPoints; ++i )
657  {
658  qreal x, y;
659  t.map( mX.at( i ), mY.at( i ), &x, &y );
660  mX[i] = x;
661  mY[i] = y;
662  }
663 }
664 
665 #if 0
666 void QgsCircularStringV2::clip( const QgsRectangle& rect )
667 {
668  //todo...
669 }
670 #endif
671 
673 {
674  int nPoints = numPoints();
675  if ( nPoints < 1 )
676  {
677  return;
678  }
679 
680  if ( path.isEmpty() || path.currentPosition() != QPointF( mX[0], mY[0] ) )
681  {
682  path.moveTo( QPointF( mX[0], mY[0] ) );
683  }
684 
685  for ( int i = 0; i < ( nPoints - 2 ) ; i += 2 )
686  {
688  segmentize( QgsPointV2( mX[i], mY[i] ), QgsPointV2( mX[i + 1], mY[i + 1] ), QgsPointV2( mX[i + 2], mY[i + 2] ), pt );
689  for ( int j = 1; j < pt.size(); ++j )
690  {
691  path.lineTo( pt.at( j ).x(), pt.at( j ).y() );
692  }
693  //arcTo( path, QPointF( mX[i], mY[i] ), QPointF( mX[i + 1], mY[i + 1] ), QPointF( mX[i + 2], mY[i + 2] ) );
694  }
695 
696  //if number of points is even, connect to last point with straight line (even though the circular string is not valid)
697  if ( nPoints % 2 == 0 )
698  {
699  path.lineTo( mX[ nPoints - 1 ], mY[ nPoints - 1 ] );
700  }
701 }
702 
703 void QgsCircularStringV2::arcTo( QPainterPath& path, QPointF pt1, QPointF pt2, QPointF pt3 )
704 {
705  double centerX, centerY, radius;
706  QgsGeometryUtils::circleCenterRadius( QgsPointV2( pt1.x(), pt1.y() ), QgsPointV2( pt2.x(), pt2.y() ), QgsPointV2( pt3.x(), pt3.y() ),
707  radius, centerX, centerY );
708 
709  double p1Angle = QgsGeometryUtils::ccwAngle( pt1.y() - centerY, pt1.x() - centerX );
710  double sweepAngle = QgsGeometryUtils::sweepAngle( centerX, centerY, pt1.x(), pt1.y(), pt2.x(), pt2.y(), pt3.x(), pt3.y() );
711 
712  double diameter = 2 * radius;
713  path.arcTo( centerX - radius, centerY - radius, diameter, diameter, p1Angle, sweepAngle );
714 }
715 
717 {
718  draw( p );
719 }
720 
722 {
723  if ( position.vertex > mX.size() || position.vertex < 1 )
724  {
725  return false;
726  }
727 
728  mX.insert( position.vertex, vertex.x() );
729  mY.insert( position.vertex, vertex.y() );
730  if ( is3D() )
731  {
732  mZ.insert( position.vertex, vertex.z() );
733  }
734  if ( isMeasure() )
735  {
736  mM.insert( position.vertex, vertex.m() );
737  }
738 
739  bool vertexNrEven = ( position.vertex % 2 == 0 );
740  if ( vertexNrEven )
741  {
742  insertVertexBetween( position.vertex - 2, position.vertex - 1, position.vertex );
743  }
744  else
745  {
746  insertVertexBetween( position.vertex, position.vertex + 1, position.vertex - 1 );
747  }
748  clearCache(); //set bounding box invalid
749  return true;
750 }
751 
753 {
754  if ( position.vertex < 0 || position.vertex >= mX.size() )
755  {
756  return false;
757  }
758 
759  mX[position.vertex] = newPos.x();
760  mY[position.vertex] = newPos.y();
761  if ( is3D() && newPos.is3D() )
762  {
763  mZ[position.vertex] = newPos.z();
764  }
765  if ( isMeasure() && newPos.isMeasure() )
766  {
767  mM[position.vertex] = newPos.m();
768  }
769  clearCache(); //set bounding box invalid
770  return true;
771 }
772 
774 {
775  int nVertices = this->numPoints();
776  if ( nVertices < 4 ) //circular string must have at least 3 vertices
777  {
778  clear();
779  return true;
780  }
781  if ( position.vertex < 0 || position.vertex > ( nVertices - 1 ) )
782  {
783  return false;
784  }
785 
786  if ( position.vertex < ( nVertices - 2 ) )
787  {
788  //remove this and the following vertex
789  deleteVertex( position.vertex + 1 );
790  deleteVertex( position.vertex );
791  }
792  else //remove this and the preceding vertex
793  {
794  deleteVertex( position.vertex );
795  deleteVertex( position.vertex - 1 );
796  }
797 
798  clearCache(); //set bounding box invalid
799  return true;
800 }
801 
803 {
804  mX.remove( i );
805  mY.remove( i );
806  if ( is3D() )
807  {
808  mZ.remove( i );
809  }
810  if ( isMeasure() )
811  {
812  mM.remove( i );
813  }
814  clearCache();
815 }
816 
817 double QgsCircularStringV2::closestSegment( const QgsPointV2& pt, QgsPointV2& segmentPt, QgsVertexId& vertexAfter, bool* leftOf, double epsilon ) const
818 {
819  Q_UNUSED( epsilon );
820  double minDist = std::numeric_limits<double>::max();
821  QgsPointV2 minDistSegmentPoint;
822  QgsVertexId minDistVertexAfter;
823  bool minDistLeftOf = false;
824 
825  double currentDist = 0.0;
826 
827  int nPoints = numPoints();
828  for ( int i = 0; i < ( nPoints - 2 ) ; i += 2 )
829  {
830  currentDist = closestPointOnArc( mX[i], mY[i], mX[i + 1], mY[i + 1], mX[i + 2], mY[i + 2], pt, segmentPt, vertexAfter, leftOf, epsilon );
831  if ( currentDist < minDist )
832  {
833  minDist = currentDist;
834  minDistSegmentPoint = segmentPt;
835  minDistVertexAfter.vertex = vertexAfter.vertex + i;
836  if ( leftOf )
837  {
838  minDistLeftOf = *leftOf;
839  }
840  }
841  }
842 
843  if ( minDist == std::numeric_limits<double>::max() )
844  return -1; // error: no segments
845 
846  segmentPt = minDistSegmentPoint;
847  vertexAfter = minDistVertexAfter;
848  vertexAfter.part = 0;
849  vertexAfter.ring = 0;
850  if ( leftOf )
851  {
852  *leftOf = minDistLeftOf;
853  }
854  return minDist;
855 }
856 
858 {
859  if ( node >= numPoints() )
860  {
861  return false;
862  }
863  point = pointN( node );
864  type = ( node % 2 == 0 ) ? QgsVertexId::SegmentVertex : QgsVertexId::CurveVertex;
865  return true;
866 }
867 
868 void QgsCircularStringV2::sumUpArea( double& sum ) const
869 {
870  int maxIndex = numPoints() - 1;
871 
872  for ( int i = 0; i < maxIndex; i += 2 )
873  {
874  QgsPointV2 p1( mX[i], mY[i] );
875  QgsPointV2 p2( mX[i + 1], mY[i + 1] );
876  QgsPointV2 p3( mX[i + 2], mY[i + 2] );
877 
878  //segment is a full circle, p2 is the center point
879  if ( p1 == p3 )
880  {
881  double r2 = QgsGeometryUtils::sqrDistance2D( p1, p2 ) / 4.0;
882  sum += M_PI * r2;
883  continue;
884  }
885 
886  sum += 0.5 * ( mX[i] * mY[i+2] - mY[i] * mX[i+2] );
887 
888  //calculate area between circle and chord, then sum / subtract from total area
889  double midPointX = ( p1.x() + p3.x() ) / 2.0;
890  double midPointY = ( p1.y() + p3.y() ) / 2.0;
891 
892  double radius, centerX, centerY;
893  QgsGeometryUtils::circleCenterRadius( p1, p2, p3, radius, centerX, centerY );
894 
895  double d = sqrt( QgsGeometryUtils::sqrDistance2D( QgsPointV2( centerX, centerY ), QgsPointV2( midPointX, midPointY ) ) );
896  double r2 = radius * radius;
897 
898  if ( d > radius )
899  {
900  //d cannot be greater than radius, something must be wrong...
901  continue;
902  }
903 
904  bool circlePointLeftOfLine = QgsGeometryUtils::leftOfLine( p2.x(), p2.y(), p1.x(), p1.y(), p3.x(), p3.y() ) < 0;
905  bool centerPointLeftOfLine = QgsGeometryUtils::leftOfLine( centerX, centerY, p1.x(), p1.y(), p3.x(), p3.y() ) < 0;
906 
907  double cov = 0.5 - d * sqrt( r2 - d * d ) / ( M_PI * r2 ) - 1 / M_PI * asin( d / radius );
908  double circleChordArea = 0;
909  if ( circlePointLeftOfLine == centerPointLeftOfLine )
910  {
911  circleChordArea = M_PI * r2 * ( 1 - cov );
912  }
913  else
914  {
915  circleChordArea = M_PI * r2 * cov;
916  }
917 
918  if ( !circlePointLeftOfLine )
919  {
920  sum += circleChordArea;
921  }
922  else
923  {
924  sum -= circleChordArea;
925  }
926  }
927 }
928 
929 double QgsCircularStringV2::closestPointOnArc( double x1, double y1, double x2, double y2, double x3, double y3,
930  const QgsPointV2& pt, QgsPointV2& segmentPt, QgsVertexId& vertexAfter, bool* leftOf, double epsilon )
931 {
932  double radius, centerX, centerY;
933  QgsPointV2 pt1( x1, y1 );
934  QgsPointV2 pt2( x2, y2 );
935  QgsPointV2 pt3( x3, y3 );
936 
937  QgsGeometryUtils::circleCenterRadius( pt1, pt2, pt3, radius, centerX, centerY );
938  double angle = QgsGeometryUtils::ccwAngle( pt.y() - centerY, pt.x() - centerX );
939  double angle1 = QgsGeometryUtils::ccwAngle( pt1.y() - centerY, pt1.x() - centerX );
940  double angle2 = QgsGeometryUtils::ccwAngle( pt2.y() - centerY, pt2.x() - centerX );
941  double angle3 = QgsGeometryUtils::ccwAngle( pt3.y() - centerY, pt3.x() - centerX );
942 
943  bool clockwise = QgsGeometryUtils::circleClockwise( angle1, angle2, angle3 );
944 
945  if ( QgsGeometryUtils::angleOnCircle( angle, angle1, angle2, angle3 ) )
946  {
947  //get point on line center -> pt with distance radius
948  segmentPt = QgsGeometryUtils::pointOnLineWithDistance( QgsPointV2( centerX, centerY ), pt, radius );
949 
950  //vertexAfter
951  vertexAfter.vertex = QgsGeometryUtils::circleAngleBetween( angle, angle1, angle2, clockwise ) ? 1 : 2;
952  }
953  else
954  {
955  double distPtPt1 = QgsGeometryUtils::sqrDistance2D( pt, pt1 );
956  double distPtPt3 = QgsGeometryUtils::sqrDistance2D( pt, pt3 );
957  segmentPt = ( distPtPt1 <= distPtPt3 ) ? pt1 : pt3;
958  vertexAfter.vertex = ( distPtPt1 <= distPtPt3 ) ? 1 : 2;
959  }
960 
961  double sqrDistance = QgsGeometryUtils::sqrDistance2D( segmentPt, pt );
962  //prevent rounding errors if the point is directly on the segment
963  if ( qgsDoubleNear( sqrDistance, 0.0, epsilon ) )
964  {
965  segmentPt.setX( pt.x() );
966  segmentPt.setY( pt.y() );
967  sqrDistance = 0.0;
968  }
969 
970  if ( leftOf )
971  {
972  *leftOf = clockwise ? sqrDistance > radius : sqrDistance < radius;
973  }
974 
975  return sqrDistance;
976 }
977 
978 void QgsCircularStringV2::insertVertexBetween( int after, int before, int pointOnCircle )
979 {
980  double xAfter = mX.at( after );
981  double yAfter = mY.at( after );
982  double xBefore = mX.at( before );
983  double yBefore = mY.at( before );
984  double xOnCircle = mX.at( pointOnCircle );
985  double yOnCircle = mY.at( pointOnCircle );
986 
987  double radius, centerX, centerY;
988  QgsGeometryUtils::circleCenterRadius( QgsPointV2( xAfter, yAfter ), QgsPointV2( xBefore, yBefore ), QgsPointV2( xOnCircle, yOnCircle ), radius, centerX, centerY );
989 
990  double x = ( xAfter + xBefore ) / 2.0;
991  double y = ( yAfter + yBefore ) / 2.0;
992 
993  QgsPointV2 newVertex = QgsGeometryUtils::pointOnLineWithDistance( QgsPointV2( centerX, centerY ), QgsPointV2( x, y ), radius );
994  mX.insert( before, newVertex.x() );
995  mY.insert( before, newVertex.y() );
996 
997  if ( is3D() )
998  {
999  mZ.insert( before, ( mZ[after] + mZ[before] ) / 2.0 );
1000  }
1001  if ( isMeasure() )
1002  {
1003  mM.insert( before, ( mM[after] + mM[before] ) / 2.0 );
1004  }
1005  clearCache();
1006 }
1007 
1009 {
1010  int before = vId.vertex - 1;
1011  int vertex = vId.vertex;
1012  int after = vId.vertex + 1;
1013 
1014  if ( vId.vertex % 2 != 0 ) // a curve vertex
1015  {
1016  if ( vId.vertex >= 1 && vId.vertex < numPoints() - 1 )
1017  {
1018  return QgsGeometryUtils::circleTangentDirection( QgsPointV2( mX[vertex], mY[vertex] ), QgsPointV2( mX[before], mY[before] ),
1019  QgsPointV2( mX[vertex], mY[vertex] ), QgsPointV2( mX[after], mY[after] ) );
1020  }
1021  }
1022  else //a point vertex
1023  {
1024  if ( vId.vertex == 0 )
1025  {
1026  return QgsGeometryUtils::circleTangentDirection( QgsPointV2( mX[0], mY[0] ), QgsPointV2( mX[0], mY[0] ),
1027  QgsPointV2( mX[1], mY[1] ), QgsPointV2( mX[2], mY[2] ) );
1028  }
1029  if ( vId.vertex >= numPoints() - 1 )
1030  {
1031  if ( numPoints() < 3 )
1032  {
1033  return 0.0;
1034  }
1035  int a = numPoints() - 3;
1036  int b = numPoints() - 2;
1037  int c = numPoints() - 1;
1038  return QgsGeometryUtils::circleTangentDirection( QgsPointV2( mX[c], mY[c] ), QgsPointV2( mX[a], mY[a] ),
1039  QgsPointV2( mX[b], mY[b] ), QgsPointV2( mX[c], mY[c] ) );
1040  }
1041  else
1042  {
1043  if ( vId.vertex + 2 > numPoints() - 1 )
1044  {
1045  return 0.0;
1046  }
1047 
1048  int vertex1 = vId.vertex - 2;
1049  int vertex2 = vId.vertex - 1;
1050  int vertex3 = vId.vertex;
1051  double angle1 = QgsGeometryUtils::circleTangentDirection( QgsPointV2( mX[vertex3], mY[vertex3] ),
1052  QgsPointV2( mX[vertex1], mY[vertex1] ), QgsPointV2( mX[vertex2], mY[vertex2] ), QgsPointV2( mX[vertex3], mY[vertex3] ) );
1053  int vertex4 = vId.vertex + 1;
1054  int vertex5 = vId.vertex + 2;
1055  double angle2 = QgsGeometryUtils::circleTangentDirection( QgsPointV2( mX[vertex3], mY[vertex3] ),
1056  QgsPointV2( mX[vertex3], mY[vertex3] ), QgsPointV2( mX[vertex4], mY[vertex4] ), QgsPointV2( mX[vertex5], mY[vertex5] ) );
1057  return QgsGeometryUtils::averageAngle( angle1, angle2 );
1058  }
1059  }
1060  return 0.0;
1061 }
1062 
1064 {
1065  QgsCircularStringV2* copy = clone();
1066  std::reverse( copy->mX.begin(), copy->mX.end() );
1067  std::reverse( copy->mY.begin(), copy->mY.end() );
1068  if ( is3D() )
1069  {
1070  std::reverse( copy->mZ.begin(), copy->mZ.end() );
1071  }
1072  if ( isMeasure() )
1073  {
1074  std::reverse( copy->mM.begin(), copy->mM.end() );
1075  }
1076  return copy;
1077 }
1078 
1079 bool QgsCircularStringV2::addZValue( double zValue )
1080 {
1081  if ( QgsWKBTypes::hasZ( mWkbType ) )
1082  return false;
1083 
1084  clearCache();
1086 
1087  int nPoints = numPoints();
1088  mZ.clear();
1089  mZ.reserve( nPoints );
1090  for ( int i = 0; i < nPoints; ++i )
1091  {
1092  mZ << zValue;
1093  }
1094  return true;
1095 }
1096 
1097 bool QgsCircularStringV2::addMValue( double mValue )
1098 {
1099  if ( QgsWKBTypes::hasM( mWkbType ) )
1100  return false;
1101 
1102  clearCache();
1104 
1105  int nPoints = numPoints();
1106  mM.clear();
1107  mM.reserve( nPoints );
1108  for ( int i = 0; i < nPoints; ++i )
1109  {
1110  mM << mValue;
1111  }
1112  return true;
1113 }
1114 
1116 {
1117  if ( !QgsWKBTypes::hasZ( mWkbType ) )
1118  return false;
1119 
1120  clearCache();
1121 
1123  mZ.clear();
1124  return true;
1125 }
1126 
1128 {
1129  if ( !QgsWKBTypes::hasM( mWkbType ) )
1130  return false;
1131 
1132  clearCache();
1133 
1135  mM.clear();
1136  return true;
1137 }
QString wktTypeStr() const
Returns the WKT type string of the geometry.
virtual bool dropMValue() override
Drops any measure values which exist in the geometry.
void clear()
void transformCoords(int numPoint, double *x, double *y, double *z, TransformDirection direction=ForwardTransform) const
Transform an array of coordinates to a different Coordinate System If the direction is ForwardTransfo...
A rectangle specified with double values.
Definition: qgsrectangle.h:35
static QgsPointV2 pointOnLineWithDistance(const QgsPointV2 &startPoint, const QgsPointV2 &directionPoint, double distance)
Returns a point a specified distance toward a second point.
QString asWkt(int precision=17) const override
Returns a WKT representation of the geometry.
QDomElement asGML2(QDomDocument &doc, int precision=17, const QString &ns="gml") const override
Returns a GML2 representation of the geometry.
static double circleTangentDirection(const QgsPointV2 &tangentPoint, const QgsPointV2 &cp1, const QgsPointV2 &cp2, const QgsPointV2 &cp3)
Calculates the direction angle of a circle tangent (clockwise from north in radians) ...
virtual QgsRectangle calculateBoundingBox() const override
Default calculator for the minimal bounding box for the geometry.
static QPair< QgsWKBTypes::Type, QString > wktReadBlock(const QString &wkt)
Parses a WKT block of the format "TYPE( contents )" and returns a pair of geometry type to contents (...
int wkbSize() const override
Returns the size of the WKB representation of the geometry.
static bool circleAngleBetween(double angle, double angle1, double angle2, bool clockwise)
Returns true if, in a circle, angle is between angle1 and angle2.
QPointF currentPosition() const
static double ccwAngle(double dy, double dx)
Returns the counter clockwise angle between a line with components dx, dy and the line with dx > 0 an...
virtual QgsLineStringV2 * curveToLine(double tolerance=M_PI_2/90, SegmentationToleranceType toleranceType=MaximumAngle) const override
Returns a new line string geometry corresponding to a segmentized approximation of the curve...
QDomNode appendChild(const QDomNode &newChild)
iterator begin()
void push_back(const T &value)
static double averageAngle(double x1, double y1, double x2, double y2, double x3, double y3)
Angle between two linear segments.
static Type addZ(Type type)
Adds the z dimension to a WKB type and returns the new type.
Definition: qgswkbtypes.h:757
static bool hasM(Type type)
Tests whether a WKB type contains m values.
Definition: qgswkbtypes.h:714
QPoint map(const QPoint &point) const
Circular string geometry type.
static void pointsToWKB(QgsWkbPtr &wkb, const QgsPointSequenceV2 &points, bool is3D, bool isMeasure)
Returns a LinearRing { uint32 numPoints; Point points[numPoints]; }.
virtual bool fromWkt(const QString &wkt) override
Sets the geometry from a WKT string.
virtual bool addMValue(double mValue=0) override
Adds a measure to the geometry, initialized to a preset value.
SegmentationToleranceType
Segmentation tolerance as maximum angle or maximum difference between approximation and circle...
TransformDirection
Enum used to indicate the direction (forward or inverse) of the transform.
const T & at(int i) const
void addToPainterPath(QPainterPath &path) const override
Adds a curve to a painter path.
void insert(int i, const T &value)
static void circleCenterRadius(const QgsPointV2 &pt1, const QgsPointV2 &pt2, const QgsPointV2 &pt3, double &radius, double &centerX, double &centerY)
Returns radius and center of the circle through pt1, pt2, pt3.
virtual bool operator!=(const QgsCurveV2 &other) const override
void moveTo(const QPointF &point)
void setX(double x)
Sets the point&#39;s x-coordinate.
Definition: qgspointv2.h:124
static double leftOfLine(double x, double y, double x1, double y1, double x2, double y2)
Returns < 0 if point(x/y) is left of the line x1,y1 -> x2,y2.
virtual void clear() override
Clears the geometry, ie reset it to a null geometry.
static double circleLength(double x1, double y1, double x2, double y2, double x3, double y3)
Length of a circular string segment defined by pt1, pt2, pt3.
void transform(const QgsCoordinateTransform &ct, QgsCoordinateTransform::TransformDirection d=QgsCoordinateTransform::ForwardTransform, bool transformZ=false) override
Transforms the geometry using a coordinate transform.
QDomElement createElementNS(const QString &nsURI, const QString &qName)
double z() const
Returns the point&#39;s z-coordinate.
Definition: qgspointv2.h:80
static bool hasZ(Type type)
Tests whether a WKB type contains the z-dimension.
Definition: qgswkbtypes.h:667
double y() const
Returns the point&#39;s y-coordinate.
Definition: qgspointv2.h:74
static endian_t endian()
Returns whether this machine uses big or little endian.
virtual bool fromWkb(QgsConstWkbPtr wkb) override
Sets the geometry from a WKB string.
bool qgsDoubleNear(double a, double b, double epsilon=4 *DBL_EPSILON)
Compare two doubles (but allow some difference)
Definition: qgis.h:353
void setY(double y)
Sets the point&#39;s y-coordinate.
Definition: qgspointv2.h:130
int size() const
virtual double length() const override
Returns the length of the geometry.
double ANALYSIS_EXPORT max(double x, double y)
Returns the maximum of two doubles or the first argument if both are equal.
T * data()
virtual bool deleteVertex(QgsVertexId position) override
Deletes a vertex within the geometry.
void clear()
static bool circleClockwise(double angle1, double angle2, double angle3)
Returns true if circle is ordered clockwise.
qreal x() const
qreal y() const
void append(const T &value)
static Type dropZ(Type type)
Drops the z dimension (if present) for a WKB type and returns the new type.
Definition: qgswkbtypes.h:811
void drawAsPolygon(QPainter &p) const override
Draws the curve as a polygon on the specified QPainter.
static Type addM(Type type)
Adds the m dimension to a WKB type and returns the new type.
Definition: qgswkbtypes.h:781
void resize(int size)
virtual void clearCache() const override
Clears any cached parameters associated with the geometry, eg bounding boxes.
Definition: qgscurvev2.h:120
static double sqrDistance2D(const QgsPointV2 &pt1, const QgsPointV2 &pt2)
Returns the squared 2D distance between two points.
Utility class for identifying a unique vertex within a geometry.
QDomElement asGML2(QDomDocument &doc, int precision=17, const QString &ns="gml") const override
Returns a GML2 representation of the geometry.
Line string geometry type, with support for z-dimension and m-values.
static QString pointsToWKT(const QgsPointSequenceV2 &points, int precision, bool is3D, bool isMeasure)
Returns a WKT coordinate list.
void lineTo(const QPointF &endPoint)
void setPoints(const QgsPointSequenceV2 &points)
Resets the line string to match the specified list of points.
bool isMeasure() const
Returns true if the geometry contains m values.
Point geometry type, with support for z-dimension and m-values.
Definition: qgspointv2.h:34
QgsPointV2 pointN(int i) const
Returns the point at index i within the circular string.
bool pointAt(int node, QgsPointV2 &point, QgsVertexId::VertexType &type) const override
Returns the point and vertex id of a point within the curve.
void remove(int i)
#define M_PI
void sumUpArea(double &sum) const override
Sums up the area of the curve by iterating over the vertices (shoelace formula).
virtual bool moveVertex(QgsVertexId position, const QgsPointV2 &newPos) override
Moves a vertex within the geometry.
double x() const
Returns the point&#39;s x-coordinate.
Definition: qgspointv2.h:68
void setZMTypeFromSubGeometry(const QgsAbstractGeometryV2 *subggeom, QgsWKBTypes::Type baseGeomType)
Updates the geometry type based on whether sub geometries contain z or m values.
static bool angleOnCircle(double angle, double angle1, double angle2, double angle3)
Returns true if an angle is between angle1 and angle3 on a circle described by angle1, angle2 and angle3.
virtual QgsPointV2 startPoint() const override
Returns the starting point of the curve.
double closestSegment(const QgsPointV2 &pt, QgsPointV2 &segmentPt, QgsVertexId &vertexAfter, bool *leftOf, double epsilon) const override
Searches for the closest segment of the geometry to a given point.
QString asJSON(int precision=17) const override
Returns a GeoJSON representation of the geometry.
virtual QgsCircularStringV2 * reversed() const override
Returns a reversed copy of the curve, where the direction of the curve has been flipped.
static QDomElement pointsToGML3(const QgsPointSequenceV2 &points, QDomDocument &doc, int precision, const QString &ns, bool is3D)
Returns a gml::posList DOM element.
unsigned char * asWkb(int &binarySize) const override
Returns a WKB representation of the geometry.
double ANALYSIS_EXPORT angle(Point3D *p1, Point3D *p2, Point3D *p3, Point3D *p4)
Calculates the angle between two segments (in 2 dimension, z-values are ignored)
void reserve(int size)
virtual bool addZValue(double zValue=0) override
Adds a z-dimension to the geometry, initialized to a preset value.
virtual QgsCircularStringV2 * clone() const override
Clones the geometry by performing a deep copy.
void combineExtentWith(const QgsRectangle &rect)
expand the rectangle so that covers both the original rectangle and the given rectangle ...
void setPoints(const QgsPointSequenceV2 &points)
Sets the circular string&#39;s points.
virtual QgsPointV2 endPoint() const override
Returns the end point of the curve.
static Type dropM(Type type)
Drops the m dimension (if present) for a WKB type and returns the new type.
Definition: qgswkbtypes.h:828
const T & at(int i) const
virtual bool dropZValue() override
Drops any z-dimensions which exist in the geometry.
void drawPath(const QPainterPath &path)
static QgsPointSequenceV2 pointsFromWKT(const QString &wktCoordinateList, bool is3D, bool isMeasure)
Returns a list of points contained in a WKT string.
static double sweepAngle(double centerX, double centerY, double x1, double y1, double x2, double y2, double x3, double y3)
Calculates angle of a circular string part defined by pt1, pt2, pt3.
void insert(int i, const T &value)
QgsWKBTypes::Type wkbType() const
Returns the WKB type of the geometry.
void draw(QPainter &p) const override
Draws the geometry using the specified QPainter.
bool isEmpty() const
Class for doing transforms between two map coordinate systems.
double vertexAngle(QgsVertexId vertex) const override
Returns approximate rotation angle for a vertex.
static Type flatType(Type type)
Returns the flat type for a WKB type.
Definition: qgswkbtypes.h:366
void points(QgsPointSequenceV2 &pts) const override
Returns a list of points within the curve.
double ANALYSIS_EXPORT leftOf(Point3D *thepoint, Point3D *p1, Point3D *p2)
Returns whether &#39;thepoint&#39; is left or right of the line from &#39;p1&#39; to &#39;p2&#39;.
QgsWKBTypes::Type readHeader() const
Definition: qgswkbptr.cpp:38
const_iterator constEnd() const
const_iterator constBegin() const
Abstract base class for curved geometry type.
Definition: qgscurvev2.h:32
QDomElement asGML3(QDomDocument &doc, int precision=17, const QString &ns="gml") const override
Returns a GML3 representation of the geometry.
int size() const
virtual bool insertVertex(QgsVertexId position, const QgsPointV2 &vertex) override
Inserts a vertex into the geometry.
void arcTo(const QRectF &rectangle, qreal startAngle, qreal sweepLength)
iterator end()
double m() const
Returns the point&#39;s m value.
Definition: qgspointv2.h:86
QString asJSON(int precision=17) const override
Returns a GeoJSON representation of the geometry.
virtual bool operator==(const QgsCurveV2 &other) const override
bool is3D() const
Returns true if the geometry is 3D and contains a z-value.
int numPoints() const override
Returns the number of points in the curve.