QGIS API Documentation  2.18.3-Las Palmas (77b8c3d)
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  //adapted code from postgis
479  double radius = 0;
480  double centerX = 0;
481  double centerY = 0;
482  QgsGeometryUtils::circleCenterRadius( p1, p2, p3, radius, centerX, centerY );
483  int segSide = segmentSide( p1, p3, p2 );
484 
485  if ( p1 != p3 && ( radius < 0 || qgsDoubleNear( segSide, 0.0 ) ) ) //points are colinear
486  {
487  points.append( p1 );
488  points.append( p2 );
489  points.append( p3 );
490  return;
491  }
492 
493  bool clockwise = false;
494  if ( segSide == -1 )
495  {
496  clockwise = true;
497  }
498 
499  double increment = tolerance; //one segment per degree
500  if ( toleranceType == QgsAbstractGeometryV2::MaximumDifference )
501  {
502  double halfAngle = acos( -tolerance / radius + 1 );
503  increment = 2 * halfAngle;
504  }
505 
506  //angles of pt1, pt2, pt3
507  double a1 = atan2( p1.y() - centerY, p1.x() - centerX );
508  double a2 = atan2( p2.y() - centerY, p2.x() - centerX );
509  double a3 = atan2( p3.y() - centerY, p3.x() - centerX );
510 
511  if ( clockwise )
512  {
513  increment *= -1;
514  /* Adjust a3 down so we can decrement from a1 to a3 cleanly */
515  if ( a3 >= a1 )
516  a3 -= 2.0 * M_PI;
517  if ( a2 > a1 )
518  a2 -= 2.0 * M_PI;
519  }
520  else
521  {
522  /* Adjust a3 up so we can increment from a1 to a3 cleanly */
523  if ( a3 <= a1 )
524  a3 += 2.0 * M_PI;
525  if ( a2 < a1 )
526  a2 += 2.0 * M_PI;
527  }
528 
529  bool hasZ = is3D();
530  bool hasM = isMeasure();
531 
532  double x, y;
533  double z = 0;
534  double m = 0;
535 
536  points.append( p1 );
537  if ( p2 != p3 && p1 != p2 ) //draw straight line segment if two points have the same position
538  {
539  QgsWKBTypes::Type pointWkbType = QgsWKBTypes::Point;
540  if ( hasZ )
541  pointWkbType = QgsWKBTypes::addZ( pointWkbType );
542  if ( hasM )
543  pointWkbType = QgsWKBTypes::addM( pointWkbType );
544 
545  //make sure the curve point p2 is part of the segmentized vertices. But only if p1 != p3
546  bool addP2 = true;
547  if ( qgsDoubleNear( p1.x(), p3.x() ) && qgsDoubleNear( p1.y(), p3.y() ) )
548  {
549  addP2 = false;
550  }
551 
552  for ( double angle = a1 + increment; clockwise ? angle > a3 : angle < a3; angle += increment )
553  {
554  if (( addP2 && clockwise && angle < a2 ) || ( addP2 && !clockwise && angle > a2 ) )
555  {
556  points.append( p2 );
557  addP2 = false;
558  }
559 
560  x = centerX + radius * cos( angle );
561  y = centerY + radius * sin( angle );
562 
563  if ( !hasZ && !hasM )
564  {
565  points.append( QgsPointV2( x, y ) );
566  continue;
567  }
568 
569  if ( hasZ )
570  {
571  z = interpolateArc( angle, a1, a2, a3, p1.z(), p2.z(), p3.z() );
572  }
573  if ( hasM )
574  {
575  m = interpolateArc( angle, a1, a2, a3, p1.m(), p2.m(), p3.m() );
576  }
577 
578  points.append( QgsPointV2( pointWkbType, x, y, z, m ) );
579  }
580  }
581  points.append( p3 );
582 }
583 
584 int QgsCircularStringV2::segmentSide( const QgsPointV2& pt1, const QgsPointV2& pt3, const QgsPointV2& pt2 ) const
585 {
586  double side = (( pt2.x() - pt1.x() ) * ( pt3.y() - pt1.y() ) - ( pt3.x() - pt1.x() ) * ( pt2.y() - pt1.y() ) );
587  if ( side == 0.0 )
588  {
589  return 0;
590  }
591  else
592  {
593  if ( side < 0 )
594  {
595  return -1;
596  }
597  if ( side > 0 )
598  {
599  return 1;
600  }
601  return 0;
602  }
603 }
604 
605 double QgsCircularStringV2::interpolateArc( double angle, double a1, double a2, double a3, double zm1, double zm2, double zm3 ) const
606 {
607  /* Counter-clockwise sweep */
608  if ( a1 < a2 )
609  {
610  if ( angle <= a2 )
611  return zm1 + ( zm2 - zm1 ) * ( angle - a1 ) / ( a2 - a1 );
612  else
613  return zm2 + ( zm3 - zm2 ) * ( angle - a2 ) / ( a3 - a2 );
614  }
615  /* Clockwise sweep */
616  else
617  {
618  if ( angle >= a2 )
619  return zm1 + ( zm2 - zm1 ) * ( a1 - angle ) / ( a1 - a2 );
620  else
621  return zm2 + ( zm3 - zm2 ) * ( a2 - angle ) / ( a2 - a3 );
622  }
623 }
624 
626 {
627  QPainterPath path;
628  addToPainterPath( path );
629  p.drawPath( path );
630 }
631 
633 {
634  clearCache();
635 
636  double* zArray = mZ.data();
637 
638  bool hasZ = is3D();
639  int nPoints = numPoints();
640  bool useDummyZ = !hasZ || !transformZ;
641  if ( useDummyZ )
642  {
643  zArray = new double[nPoints];
644  for ( int i = 0; i < nPoints; ++i )
645  {
646  zArray[i] = 0;
647  }
648  }
649  ct.transformCoords( nPoints, mX.data(), mY.data(), zArray, d );
650  if ( useDummyZ )
651  {
652  delete[] zArray;
653  }
654 }
655 
657 {
658  clearCache();
659 
660  int nPoints = numPoints();
661  for ( int i = 0; i < nPoints; ++i )
662  {
663  qreal x, y;
664  t.map( mX.at( i ), mY.at( i ), &x, &y );
665  mX[i] = x;
666  mY[i] = y;
667  }
668 }
669 
670 #if 0
671 void QgsCircularStringV2::clip( const QgsRectangle& rect )
672 {
673  //todo...
674 }
675 #endif
676 
678 {
679  int nPoints = numPoints();
680  if ( nPoints < 1 )
681  {
682  return;
683  }
684 
685  if ( path.isEmpty() || path.currentPosition() != QPointF( mX[0], mY[0] ) )
686  {
687  path.moveTo( QPointF( mX[0], mY[0] ) );
688  }
689 
690  for ( int i = 0; i < ( nPoints - 2 ) ; i += 2 )
691  {
693  segmentize( QgsPointV2( mX[i], mY[i] ), QgsPointV2( mX[i + 1], mY[i + 1] ), QgsPointV2( mX[i + 2], mY[i + 2] ), pt );
694  for ( int j = 1; j < pt.size(); ++j )
695  {
696  path.lineTo( pt.at( j ).x(), pt.at( j ).y() );
697  }
698  //arcTo( path, QPointF( mX[i], mY[i] ), QPointF( mX[i + 1], mY[i + 1] ), QPointF( mX[i + 2], mY[i + 2] ) );
699  }
700 
701  //if number of points is even, connect to last point with straight line (even though the circular string is not valid)
702  if ( nPoints % 2 == 0 )
703  {
704  path.lineTo( mX[ nPoints - 1 ], mY[ nPoints - 1 ] );
705  }
706 }
707 
708 void QgsCircularStringV2::arcTo( QPainterPath& path, QPointF pt1, QPointF pt2, QPointF pt3 )
709 {
710  double centerX, centerY, radius;
711  QgsGeometryUtils::circleCenterRadius( QgsPointV2( pt1.x(), pt1.y() ), QgsPointV2( pt2.x(), pt2.y() ), QgsPointV2( pt3.x(), pt3.y() ),
712  radius, centerX, centerY );
713 
714  double p1Angle = QgsGeometryUtils::ccwAngle( pt1.y() - centerY, pt1.x() - centerX );
715  double sweepAngle = QgsGeometryUtils::sweepAngle( centerX, centerY, pt1.x(), pt1.y(), pt2.x(), pt2.y(), pt3.x(), pt3.y() );
716 
717  double diameter = 2 * radius;
718  path.arcTo( centerX - radius, centerY - radius, diameter, diameter, p1Angle, sweepAngle );
719 }
720 
722 {
723  draw( p );
724 }
725 
727 {
728  if ( position.vertex > mX.size() || position.vertex < 1 )
729  {
730  return false;
731  }
732 
733  mX.insert( position.vertex, vertex.x() );
734  mY.insert( position.vertex, vertex.y() );
735  if ( is3D() )
736  {
737  mZ.insert( position.vertex, vertex.z() );
738  }
739  if ( isMeasure() )
740  {
741  mM.insert( position.vertex, vertex.m() );
742  }
743 
744  bool vertexNrEven = ( position.vertex % 2 == 0 );
745  if ( vertexNrEven )
746  {
747  insertVertexBetween( position.vertex - 2, position.vertex - 1, position.vertex );
748  }
749  else
750  {
751  insertVertexBetween( position.vertex, position.vertex + 1, position.vertex - 1 );
752  }
753  clearCache(); //set bounding box invalid
754  return true;
755 }
756 
758 {
759  if ( position.vertex < 0 || position.vertex >= mX.size() )
760  {
761  return false;
762  }
763 
764  mX[position.vertex] = newPos.x();
765  mY[position.vertex] = newPos.y();
766  if ( is3D() && newPos.is3D() )
767  {
768  mZ[position.vertex] = newPos.z();
769  }
770  if ( isMeasure() && newPos.isMeasure() )
771  {
772  mM[position.vertex] = newPos.m();
773  }
774  clearCache(); //set bounding box invalid
775  return true;
776 }
777 
779 {
780  int nVertices = this->numPoints();
781  if ( nVertices < 4 ) //circular string must have at least 3 vertices
782  {
783  clear();
784  return true;
785  }
786  if ( position.vertex < 0 || position.vertex > ( nVertices - 1 ) )
787  {
788  return false;
789  }
790 
791  if ( position.vertex < ( nVertices - 2 ) )
792  {
793  //remove this and the following vertex
794  deleteVertex( position.vertex + 1 );
795  deleteVertex( position.vertex );
796  }
797  else //remove this and the preceding vertex
798  {
799  deleteVertex( position.vertex );
800  deleteVertex( position.vertex - 1 );
801  }
802 
803  clearCache(); //set bounding box invalid
804  return true;
805 }
806 
808 {
809  mX.remove( i );
810  mY.remove( i );
811  if ( is3D() )
812  {
813  mZ.remove( i );
814  }
815  if ( isMeasure() )
816  {
817  mM.remove( i );
818  }
819  clearCache();
820 }
821 
822 double QgsCircularStringV2::closestSegment( const QgsPointV2& pt, QgsPointV2& segmentPt, QgsVertexId& vertexAfter, bool* leftOf, double epsilon ) const
823 {
824  Q_UNUSED( epsilon );
825  double minDist = std::numeric_limits<double>::max();
826  QgsPointV2 minDistSegmentPoint;
827  QgsVertexId minDistVertexAfter;
828  bool minDistLeftOf = false;
829 
830  double currentDist = 0.0;
831 
832  int nPoints = numPoints();
833  for ( int i = 0; i < ( nPoints - 2 ) ; i += 2 )
834  {
835  currentDist = closestPointOnArc( mX[i], mY[i], mX[i + 1], mY[i + 1], mX[i + 2], mY[i + 2], pt, segmentPt, vertexAfter, leftOf, epsilon );
836  if ( currentDist < minDist )
837  {
838  minDist = currentDist;
839  minDistSegmentPoint = segmentPt;
840  minDistVertexAfter.vertex = vertexAfter.vertex + i;
841  if ( leftOf )
842  {
843  minDistLeftOf = *leftOf;
844  }
845  }
846  }
847 
848  if ( minDist == std::numeric_limits<double>::max() )
849  return -1; // error: no segments
850 
851  segmentPt = minDistSegmentPoint;
852  vertexAfter = minDistVertexAfter;
853  vertexAfter.part = 0;
854  vertexAfter.ring = 0;
855  if ( leftOf )
856  {
857  *leftOf = minDistLeftOf;
858  }
859  return minDist;
860 }
861 
863 {
864  if ( node >= numPoints() )
865  {
866  return false;
867  }
868  point = pointN( node );
869  type = ( node % 2 == 0 ) ? QgsVertexId::SegmentVertex : QgsVertexId::CurveVertex;
870  return true;
871 }
872 
873 void QgsCircularStringV2::sumUpArea( double& sum ) const
874 {
875  int maxIndex = numPoints() - 1;
876 
877  for ( int i = 0; i < maxIndex; i += 2 )
878  {
879  QgsPointV2 p1( mX[i], mY[i] );
880  QgsPointV2 p2( mX[i + 1], mY[i + 1] );
881  QgsPointV2 p3( mX[i + 2], mY[i + 2] );
882 
883  //segment is a full circle, p2 is the center point
884  if ( p1 == p3 )
885  {
886  double r2 = QgsGeometryUtils::sqrDistance2D( p1, p2 ) / 4.0;
887  sum += M_PI * r2;
888  continue;
889  }
890 
891  sum += 0.5 * ( mX[i] * mY[i+2] - mY[i] * mX[i+2] );
892 
893  //calculate area between circle and chord, then sum / subtract from total area
894  double midPointX = ( p1.x() + p3.x() ) / 2.0;
895  double midPointY = ( p1.y() + p3.y() ) / 2.0;
896 
897  double radius, centerX, centerY;
898  QgsGeometryUtils::circleCenterRadius( p1, p2, p3, radius, centerX, centerY );
899 
900  double d = sqrt( QgsGeometryUtils::sqrDistance2D( QgsPointV2( centerX, centerY ), QgsPointV2( midPointX, midPointY ) ) );
901  double r2 = radius * radius;
902 
903  if ( d > radius )
904  {
905  //d cannot be greater than radius, something must be wrong...
906  continue;
907  }
908 
909  bool circlePointLeftOfLine = QgsGeometryUtils::leftOfLine( p2.x(), p2.y(), p1.x(), p1.y(), p3.x(), p3.y() ) < 0;
910  bool centerPointLeftOfLine = QgsGeometryUtils::leftOfLine( centerX, centerY, p1.x(), p1.y(), p3.x(), p3.y() ) < 0;
911 
912  double cov = 0.5 - d * sqrt( r2 - d * d ) / ( M_PI * r2 ) - 1 / M_PI * asin( d / radius );
913  double circleChordArea = 0;
914  if ( circlePointLeftOfLine == centerPointLeftOfLine )
915  {
916  circleChordArea = M_PI * r2 * ( 1 - cov );
917  }
918  else
919  {
920  circleChordArea = M_PI * r2 * cov;
921  }
922 
923  if ( !circlePointLeftOfLine )
924  {
925  sum += circleChordArea;
926  }
927  else
928  {
929  sum -= circleChordArea;
930  }
931  }
932 }
933 
934 double QgsCircularStringV2::closestPointOnArc( double x1, double y1, double x2, double y2, double x3, double y3,
935  const QgsPointV2& pt, QgsPointV2& segmentPt, QgsVertexId& vertexAfter, bool* leftOf, double epsilon )
936 {
937  double radius, centerX, centerY;
938  QgsPointV2 pt1( x1, y1 );
939  QgsPointV2 pt2( x2, y2 );
940  QgsPointV2 pt3( x3, y3 );
941 
942  QgsGeometryUtils::circleCenterRadius( pt1, pt2, pt3, radius, centerX, centerY );
943  double angle = QgsGeometryUtils::ccwAngle( pt.y() - centerY, pt.x() - centerX );
944  double angle1 = QgsGeometryUtils::ccwAngle( pt1.y() - centerY, pt1.x() - centerX );
945  double angle2 = QgsGeometryUtils::ccwAngle( pt2.y() - centerY, pt2.x() - centerX );
946  double angle3 = QgsGeometryUtils::ccwAngle( pt3.y() - centerY, pt3.x() - centerX );
947 
948  bool clockwise = QgsGeometryUtils::circleClockwise( angle1, angle2, angle3 );
949 
950  if ( QgsGeometryUtils::angleOnCircle( angle, angle1, angle2, angle3 ) )
951  {
952  //get point on line center -> pt with distance radius
953  segmentPt = QgsGeometryUtils::pointOnLineWithDistance( QgsPointV2( centerX, centerY ), pt, radius );
954 
955  //vertexAfter
956  vertexAfter.vertex = QgsGeometryUtils::circleAngleBetween( angle, angle1, angle2, clockwise ) ? 1 : 2;
957  }
958  else
959  {
960  double distPtPt1 = QgsGeometryUtils::sqrDistance2D( pt, pt1 );
961  double distPtPt3 = QgsGeometryUtils::sqrDistance2D( pt, pt3 );
962  segmentPt = ( distPtPt1 <= distPtPt3 ) ? pt1 : pt3;
963  vertexAfter.vertex = ( distPtPt1 <= distPtPt3 ) ? 1 : 2;
964  }
965 
966  double sqrDistance = QgsGeometryUtils::sqrDistance2D( segmentPt, pt );
967  //prevent rounding errors if the point is directly on the segment
968  if ( qgsDoubleNear( sqrDistance, 0.0, epsilon ) )
969  {
970  segmentPt.setX( pt.x() );
971  segmentPt.setY( pt.y() );
972  sqrDistance = 0.0;
973  }
974 
975  if ( leftOf )
976  {
977  *leftOf = clockwise ? sqrDistance > radius : sqrDistance < radius;
978  }
979 
980  return sqrDistance;
981 }
982 
983 void QgsCircularStringV2::insertVertexBetween( int after, int before, int pointOnCircle )
984 {
985  double xAfter = mX.at( after );
986  double yAfter = mY.at( after );
987  double xBefore = mX.at( before );
988  double yBefore = mY.at( before );
989  double xOnCircle = mX.at( pointOnCircle );
990  double yOnCircle = mY.at( pointOnCircle );
991 
992  double radius, centerX, centerY;
993  QgsGeometryUtils::circleCenterRadius( QgsPointV2( xAfter, yAfter ), QgsPointV2( xBefore, yBefore ), QgsPointV2( xOnCircle, yOnCircle ), radius, centerX, centerY );
994 
995  double x = ( xAfter + xBefore ) / 2.0;
996  double y = ( yAfter + yBefore ) / 2.0;
997 
998  QgsPointV2 newVertex = QgsGeometryUtils::pointOnLineWithDistance( QgsPointV2( centerX, centerY ), QgsPointV2( x, y ), radius );
999  mX.insert( before, newVertex.x() );
1000  mY.insert( before, newVertex.y() );
1001 
1002  if ( is3D() )
1003  {
1004  mZ.insert( before, ( mZ[after] + mZ[before] ) / 2.0 );
1005  }
1006  if ( isMeasure() )
1007  {
1008  mM.insert( before, ( mM[after] + mM[before] ) / 2.0 );
1009  }
1010  clearCache();
1011 }
1012 
1014 {
1015  int before = vId.vertex - 1;
1016  int vertex = vId.vertex;
1017  int after = vId.vertex + 1;
1018 
1019  if ( vId.vertex % 2 != 0 ) // a curve vertex
1020  {
1021  if ( vId.vertex >= 1 && vId.vertex < numPoints() - 1 )
1022  {
1023  return QgsGeometryUtils::circleTangentDirection( QgsPointV2( mX[vertex], mY[vertex] ), QgsPointV2( mX[before], mY[before] ),
1024  QgsPointV2( mX[vertex], mY[vertex] ), QgsPointV2( mX[after], mY[after] ) );
1025  }
1026  }
1027  else //a point vertex
1028  {
1029  if ( vId.vertex == 0 )
1030  {
1031  return QgsGeometryUtils::circleTangentDirection( QgsPointV2( mX[0], mY[0] ), QgsPointV2( mX[0], mY[0] ),
1032  QgsPointV2( mX[1], mY[1] ), QgsPointV2( mX[2], mY[2] ) );
1033  }
1034  if ( vId.vertex >= numPoints() - 1 )
1035  {
1036  if ( numPoints() < 3 )
1037  {
1038  return 0.0;
1039  }
1040  int a = numPoints() - 3;
1041  int b = numPoints() - 2;
1042  int c = numPoints() - 1;
1043  return QgsGeometryUtils::circleTangentDirection( QgsPointV2( mX[c], mY[c] ), QgsPointV2( mX[a], mY[a] ),
1044  QgsPointV2( mX[b], mY[b] ), QgsPointV2( mX[c], mY[c] ) );
1045  }
1046  else
1047  {
1048  if ( vId.vertex + 2 > numPoints() - 1 )
1049  {
1050  return 0.0;
1051  }
1052 
1053  int vertex1 = vId.vertex - 2;
1054  int vertex2 = vId.vertex - 1;
1055  int vertex3 = vId.vertex;
1056  double angle1 = QgsGeometryUtils::circleTangentDirection( QgsPointV2( mX[vertex3], mY[vertex3] ),
1057  QgsPointV2( mX[vertex1], mY[vertex1] ), QgsPointV2( mX[vertex2], mY[vertex2] ), QgsPointV2( mX[vertex3], mY[vertex3] ) );
1058  int vertex4 = vId.vertex + 1;
1059  int vertex5 = vId.vertex + 2;
1060  double angle2 = QgsGeometryUtils::circleTangentDirection( QgsPointV2( mX[vertex3], mY[vertex3] ),
1061  QgsPointV2( mX[vertex3], mY[vertex3] ), QgsPointV2( mX[vertex4], mY[vertex4] ), QgsPointV2( mX[vertex5], mY[vertex5] ) );
1062  return QgsGeometryUtils::averageAngle( angle1, angle2 );
1063  }
1064  }
1065  return 0.0;
1066 }
1067 
1069 {
1070  QgsCircularStringV2* copy = clone();
1071  std::reverse( copy->mX.begin(), copy->mX.end() );
1072  std::reverse( copy->mY.begin(), copy->mY.end() );
1073  if ( is3D() )
1074  {
1075  std::reverse( copy->mZ.begin(), copy->mZ.end() );
1076  }
1077  if ( isMeasure() )
1078  {
1079  std::reverse( copy->mM.begin(), copy->mM.end() );
1080  }
1081  return copy;
1082 }
1083 
1084 bool QgsCircularStringV2::addZValue( double zValue )
1085 {
1086  if ( QgsWKBTypes::hasZ( mWkbType ) )
1087  return false;
1088 
1089  clearCache();
1091 
1092  int nPoints = numPoints();
1093  mZ.clear();
1094  mZ.reserve( nPoints );
1095  for ( int i = 0; i < nPoints; ++i )
1096  {
1097  mZ << zValue;
1098  }
1099  return true;
1100 }
1101 
1102 bool QgsCircularStringV2::addMValue( double mValue )
1103 {
1104  if ( QgsWKBTypes::hasM( mWkbType ) )
1105  return false;
1106 
1107  clearCache();
1109 
1110  int nPoints = numPoints();
1111  mM.clear();
1112  mM.reserve( nPoints );
1113  for ( int i = 0; i < nPoints; ++i )
1114  {
1115  mM << mValue;
1116  }
1117  return true;
1118 }
1119 
1121 {
1122  if ( !QgsWKBTypes::hasZ( mWkbType ) )
1123  return false;
1124 
1125  clearCache();
1126 
1128  mZ.clear();
1129  return true;
1130 }
1131 
1133 {
1134  if ( !QgsWKBTypes::hasM( mWkbType ) )
1135  return false;
1136 
1137  clearCache();
1138 
1140  mM.clear();
1141  return true;
1142 }
virtual bool dropMValue() override
Drops any measure values which exist in the geometry.
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...
void clear()
QgsWKBTypes::Type wkbType() const
Returns the WKB type of the geometry.
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)
double x() const
Returns the point&#39;s x-coordinate.
Definition: qgspointv2.h:68
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)
static bool hasZ(Type type)
Tests whether a WKB type contains the z-dimension.
Definition: qgswkbtypes.h:667
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.
QString wktTypeStr() const
Returns the WKT type string of the geometry.
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.
QgsPointV2 pointN(int i) const
Returns the point at index i within the circular string.
double y() const
Returns the point&#39;s y-coordinate.
Definition: qgspointv2.h:74
double ANALYSIS_EXPORT max(double x, double y)
Returns the maximum of two doubles or the first argument if both are equal.
QgsWKBTypes::Type readHeader() const
Definition: qgswkbptr.cpp:38
T * data()
virtual bool deleteVertex(QgsVertexId position) override
Deletes a vertex within the geometry.
void clear()
bool is3D() const
Returns true if the geometry is 3D and contains a z-value.
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:121
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.
bool isMeasure() const
Returns true if the geometry contains m values.
double z() const
Returns the point&#39;s z-coordinate.
Definition: qgspointv2.h:80
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.
Point geometry type, with support for z-dimension and m-values.
Definition: qgspointv2.h:34
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
Calculates the area of the curve.
virtual bool moveVertex(QgsVertexId position, const QgsPointV2 &newPos) override
Moves a vertex within the geometry.
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 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.
double m() const
Returns the point&#39;s m value.
Definition: qgspointv2.h:86
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;.
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()
QString asJSON(int precision=17) const override
Returns a GeoJSON representation of the geometry.
virtual bool operator==(const QgsCurveV2 &other) const override
int numPoints() const override
Returns the number of points in the curve.