QGIS API Documentation  3.17.0-Master (a035f434f4)
qgscircle.cpp
Go to the documentation of this file.
1 /***************************************************************************
2  qgscircle.cpp
3  --------------
4  begin : March 2017
5  copyright : (C) 2017 by Loîc Bartoletti
6  email : lbartoletti at tuxfamily dot org
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 "qgscircle.h"
19 #include "qgslinestring.h"
20 #include "qgsgeometryutils.h"
21 #include "qgstriangle.h"
22 
23 #include <memory>
24 
26  QgsEllipse( QgsPoint(), 0.0, 0.0, 0.0 )
27 {
28 
29 }
30 
31 QgsCircle::QgsCircle( const QgsPoint &center, double radius, double azimuth ) :
32  QgsEllipse( center, radius, radius, azimuth )
33 {
34 
35 }
36 
38 {
40  double azimuth = QgsGeometryUtils::lineAngle( pt1.x(), pt1.y(), pt2.x(), pt2.y() ) * 180.0 / M_PI;
41  double radius = pt1.distance( pt2 ) / 2.0;
42 
44 
45  return QgsCircle( center, radius, azimuth );
46 }
47 
48 static bool isPerpendicular( const QgsPoint &pt1, const QgsPoint &pt2, const QgsPoint &pt3, double epsilon )
49 {
50  // check the given point are perpendicular to x or y axis
51 
52  double yDelta_a = pt2.y() - pt1.y();
53  double xDelta_a = pt2.x() - pt1.x();
54  double yDelta_b = pt3.y() - pt2.y();
55  double xDelta_b = pt3.x() - pt2.x();
56 
57  if ( ( std::fabs( xDelta_a ) <= epsilon ) && ( std::fabs( yDelta_b ) <= epsilon ) )
58  {
59  return false;
60  }
61 
62  if ( std::fabs( yDelta_a ) <= epsilon )
63  {
64  return true;
65  }
66  else if ( std::fabs( yDelta_b ) <= epsilon )
67  {
68  return true;
69  }
70  else if ( std::fabs( xDelta_a ) <= epsilon )
71  {
72  return true;
73  }
74  else if ( std::fabs( xDelta_b ) <= epsilon )
75  {
76  return true;
77  }
78 
79  return false;
80 
81 }
82 
83 QgsCircle QgsCircle::from3Points( const QgsPoint &pt1, const QgsPoint &pt2, const QgsPoint &pt3, double epsilon )
84 {
85  QgsPoint p1, p2, p3;
86 
87  if ( !isPerpendicular( pt1, pt2, pt3, epsilon ) )
88  {
89  p1 = pt1;
90  p2 = pt2;
91  p3 = pt3;
92  }
93  else if ( !isPerpendicular( pt1, pt3, pt2, epsilon ) )
94  {
95  p1 = pt1;
96  p2 = pt3;
97  p3 = pt2;
98  }
99  else if ( !isPerpendicular( pt2, pt1, pt3, epsilon ) )
100  {
101  p1 = pt2;
102  p2 = pt1;
103  p3 = pt3;
104  }
105  else if ( !isPerpendicular( pt2, pt3, pt1, epsilon ) )
106  {
107  p1 = pt2;
108  p2 = pt3;
109  p3 = pt1;
110  }
111  else if ( !isPerpendicular( pt3, pt2, pt1, epsilon ) )
112  {
113  p1 = pt3;
114  p2 = pt2;
115  p3 = pt1;
116  }
117  else if ( !isPerpendicular( pt3, pt1, pt2, epsilon ) )
118  {
119  p1 = pt3;
120  p2 = pt1;
121  p3 = pt2;
122  }
123  else
124  {
125  return QgsCircle();
126  }
128  double radius = -0.0;
129  // Paul Bourke's algorithm
130  double yDelta_a = p2.y() - p1.y();
131  double xDelta_a = p2.x() - p1.x();
132  double yDelta_b = p3.y() - p2.y();
133  double xDelta_b = p3.x() - p2.x();
134 
135  if ( qgsDoubleNear( xDelta_a, 0.0, epsilon ) || qgsDoubleNear( xDelta_b, 0.0, epsilon ) )
136  {
137  return QgsCircle();
138  }
139 
140  double aSlope = yDelta_a / xDelta_a;
141  double bSlope = yDelta_b / xDelta_b;
142 
143  // set z coordinate for center
144  QgsGeometryUtils::setZValueFromPoints( QgsPointSequence() << p1 << p2 << p3, center );
145 
146  if ( ( std::fabs( xDelta_a ) <= epsilon ) && ( std::fabs( yDelta_b ) <= epsilon ) )
147  {
148  center.setX( 0.5 * ( p2.x() + p3.x() ) );
149  center.setY( 0.5 * ( p1.y() + p2.y() ) );
150  radius = center.distance( pt1 );
151 
152  return QgsCircle( center, radius );
153  }
154 
155  if ( std::fabs( aSlope - bSlope ) <= epsilon )
156  {
157  return QgsCircle();
158  }
159 
160  center.setX(
161  ( aSlope * bSlope * ( p1.y() - p3.y() ) +
162  bSlope * ( p1.x() + p2.x() ) -
163  aSlope * ( p2.x() + p3.x() ) ) /
164  ( 2.0 * ( bSlope - aSlope ) )
165  );
166  center.setY(
167  -1.0 * ( center.x() - ( p1.x() + p2.x() ) / 2.0 ) /
168  aSlope + ( p1.y() + p2.y() ) / 2.0
169  );
170 
171  radius = center.distance( p1 );
172 
173  return QgsCircle( center, radius );
174 }
175 
177 {
178  return QgsCircle( center, diameter / 2.0, azimuth );
179 }
180 
182 {
183  double azimuth = QgsGeometryUtils::lineAngle( center.x(), center.y(), pt1.x(), pt1.y() ) * 180.0 / M_PI;
184 
185  QgsPoint centerPt( center );
186  QgsGeometryUtils::setZValueFromPoints( QgsPointSequence() << center << pt1, centerPt );
187 
188  return QgsCircle( centerPt, centerPt.distance( pt1 ), azimuth );
189 }
190 
191 QgsCircle QgsCircle::from3Tangents( const QgsPoint &pt1_tg1, const QgsPoint &pt2_tg1, const QgsPoint &pt1_tg2, const QgsPoint &pt2_tg2, const QgsPoint &pt1_tg3, const QgsPoint &pt2_tg3, double epsilon )
192 {
193  QgsPoint p1, p2, p3;
194  bool isIntersect = false;
195  QgsGeometryUtils::segmentIntersection( pt1_tg1, pt2_tg1, pt1_tg2, pt2_tg2, p1, isIntersect, epsilon );
196  if ( !isIntersect )
197  return QgsCircle();
198  QgsGeometryUtils::segmentIntersection( pt1_tg1, pt2_tg1, pt1_tg3, pt2_tg3, p2, isIntersect, epsilon );
199  if ( !isIntersect )
200  return QgsCircle();
201  QgsGeometryUtils::segmentIntersection( pt1_tg2, pt2_tg2, pt1_tg3, pt2_tg3, p3, isIntersect, epsilon );
202  if ( !isIntersect )
203  return QgsCircle();
204 
205  if ( p1.is3D() )
206  {
207  p1.convertTo( QgsWkbTypes::dropZ( p1.wkbType() ) );
208  }
209 
210  if ( p2.is3D() )
211  {
212  p2.convertTo( QgsWkbTypes::dropZ( p2.wkbType() ) );
213  }
214 
215  if ( p3.is3D() )
216  {
217  p3.convertTo( QgsWkbTypes::dropZ( p3.wkbType() ) );
218  }
219 
220  return QgsTriangle( p1, p2, p3 ).inscribedCircle();
221 }
222 
223 QgsCircle QgsCircle::minimalCircleFrom3Points( const QgsPoint &pt1, const QgsPoint &pt2, const QgsPoint &pt3, double epsilon )
224 {
225  double l1 = pt2.distance( pt3 );
226  double l2 = pt3.distance( pt1 );
227  double l3 = pt1.distance( pt2 );
228 
229  if ( ( l1 * l1 ) - ( l2 * l2 + l3 * l3 ) >= epsilon )
230  return QgsCircle().from2Points( pt2, pt3 );
231  else if ( ( l2 * l2 ) - ( l1 * l1 + l3 * l3 ) >= epsilon )
232  return QgsCircle().from2Points( pt3, pt1 );
233  else if ( ( l3 * l3 ) - ( l1 * l1 + l2 * l2 ) >= epsilon )
234  return QgsCircle().from2Points( pt1, pt2 );
235  else
236  return QgsCircle().from3Points( pt1, pt2, pt3, epsilon );
237 }
238 
239 int QgsCircle::intersections( const QgsCircle &other, QgsPoint &intersection1, QgsPoint &intersection2, bool useZ ) const
240 {
241  if ( useZ && mCenter.is3D() && other.center().is3D() && !qgsDoubleNear( mCenter.z(), other.center().z() ) )
242  return 0;
243 
244  QgsPointXY int1, int2;
245 
247  QgsPointXY( other.center() ), other.radius(),
248  int1, int2 );
249  if ( res == 0 )
250  return 0;
251 
252  intersection1 = QgsPoint( int1 );
253  intersection2 = QgsPoint( int2 );
254  if ( useZ && mCenter.is3D() )
255  {
256  intersection1.addZValue( mCenter.z() );
257  intersection2.addZValue( mCenter.z() );
258  }
259  return res;
260 }
261 
262 bool QgsCircle::tangentToPoint( const QgsPointXY &p, QgsPointXY &pt1, QgsPointXY &pt2 ) const
263 {
265 }
266 
267 int QgsCircle::outerTangents( const QgsCircle &other, QgsPointXY &line1P1, QgsPointXY &line1P2, QgsPointXY &line2P1, QgsPointXY &line2P2 ) const
268 {
270  QgsPointXY( other.center() ), other.radius(), line1P1, line1P2, line2P1, line2P2 );
271 }
272 
273 int QgsCircle::innerTangents( const QgsCircle &other, QgsPointXY &line1P1, QgsPointXY &line1P2, QgsPointXY &line2P1, QgsPointXY &line2P2 ) const
274 {
276  QgsPointXY( other.center() ), other.radius(), line1P1, line1P2, line2P1, line2P2 );
277 }
278 
280 {
281  double delta_x = std::fabs( pt1.x() - pt2.x() );
282  double delta_y = std::fabs( pt1.x() - pt2.y() );
283  if ( !qgsDoubleNear( delta_x, delta_y ) )
284  {
285  return QgsCircle();
286  }
287 
289  QgsGeometryUtils::setZValueFromPoints( QgsPointSequence() << pt1 << pt2, center );
290 
291  return QgsCircle( center, delta_x / 2.0, 0 );
292 }
293 
294 double QgsCircle::area() const
295 {
296  return M_PI * mSemiMajorAxis * mSemiMajorAxis;
297 }
298 
299 double QgsCircle::perimeter() const
300 {
301  return 2.0 * M_PI * mSemiMajorAxis;
302 }
303 
305 {
306  mSemiMajorAxis = std::fabs( semiMajorAxis );
308 }
309 
311 {
312  mSemiMajorAxis = std::fabs( semiMinorAxis );
314 }
315 
316 QVector<QgsPoint> QgsCircle::northQuadrant() const
317 {
318  QVector<QgsPoint> quad;
319  quad.append( QgsPoint( mCenter.x(), mCenter.y() + mSemiMajorAxis, mCenter.z() ) );
320  quad.append( QgsPoint( mCenter.x() + mSemiMajorAxis, mCenter.y(), mCenter.z() ) );
321  quad.append( QgsPoint( mCenter.x(), mCenter.y() - mSemiMajorAxis, mCenter.z() ) );
322  quad.append( QgsPoint( mCenter.x() - mSemiMajorAxis, mCenter.y(), mCenter.z() ) );
323 
324  return quad;
325 }
326 
328 {
329  std::unique_ptr<QgsCircularString> circString( new QgsCircularString() );
331  QVector<QgsPoint> quad;
332  if ( oriented )
333  {
334  quad = quadrant();
335  }
336  else
337  {
338  quad = northQuadrant();
339  }
340  quad.append( quad.at( 0 ) );
341  for ( QVector<QgsPoint>::const_iterator it = quad.constBegin(); it != quad.constEnd(); ++it )
342  {
343  points.append( *it );
344  }
345  circString->setPoints( points );
346 
347  return circString.release();
348 }
349 
350 bool QgsCircle::contains( const QgsPoint &point, double epsilon ) const
351 {
352  return ( mCenter.distance( point ) <= mSemiMajorAxis + epsilon );
353 }
354 
356 {
358 }
359 
360 QString QgsCircle::toString( int pointPrecision, int radiusPrecision, int azimuthPrecision ) const
361 {
362  QString rep;
363  if ( isEmpty() )
364  rep = QStringLiteral( "Empty" );
365  else
366  rep = QStringLiteral( "Circle (Center: %1, Radius: %2, Azimuth: %3)" )
367  .arg( mCenter.asWkt( pointPrecision ), 0, 's' )
368  .arg( qgsDoubleToString( mSemiMajorAxis, radiusPrecision ), 0, 'f' )
369  .arg( qgsDoubleToString( mAzimuth, azimuthPrecision ), 0, 'f' );
370 
371  return rep;
372 
373 }
double radius() const SIP_HOLDGIL
Returns the radius of the circle.
Definition: qgscircle.h:247
Circle geometry type.
Definition: qgscircle.h:43
static QgsCircle fromExtent(const QgsPoint &pt1, const QgsPoint &pt2) SIP_HOLDGIL
Constructs a circle by an extent (aka bounding box / QgsRectangle).
Definition: qgscircle.cpp:279
A rectangle specified with double values.
Definition: qgsrectangle.h:41
double y
Definition: qgspoint.h:42
static int circleCircleOuterTangents(const QgsPointXY &center1, double radius1, const QgsPointXY &center2, double radius2, QgsPointXY &line1P1, QgsPointXY &line1P2, QgsPointXY &line2P1, QgsPointXY &line2P2) SIP_HOLDGIL
Calculates the outer tangent points for two circles, centered at center1 and center2 and with radii o...
double semiMinorAxis() const SIP_HOLDGIL
Returns the semi-minor axis.
Definition: qgsellipse.h:133
static QgsCircle from3Tangents(const QgsPoint &pt1_tg1, const QgsPoint &pt2_tg1, const QgsPoint &pt1_tg2, const QgsPoint &pt2_tg2, const QgsPoint &pt1_tg3, const QgsPoint &pt2_tg3, double epsilon=1E-8) SIP_HOLDGIL
Constructs a circle by 3 tangents on the circle (aka inscribed circle of a triangle).
Definition: qgscircle.cpp:191
bool addZValue(double zValue=0) override
Adds a z-dimension to the geometry, initialized to a preset value.
Definition: qgspoint.cpp:541
void setY(double y) SIP_HOLDGIL
Sets the point&#39;s y-coordinate.
Definition: qgspoint.h:280
static int circleCircleInnerTangents(const QgsPointXY &center1, double radius1, const QgsPointXY &center2, double radius2, QgsPointXY &line1P1, QgsPointXY &line1P2, QgsPointXY &line2P1, QgsPointXY &line2P2) SIP_HOLDGIL
Calculates the inner tangent points for two circles, centered at center1 and center2 and with radii o...
static bool segmentIntersection(const QgsPoint &p1, const QgsPoint &p2, const QgsPoint &q1, const QgsPoint &q2, QgsPoint &intersectionPoint, bool &isIntersection, double tolerance=1e-8, bool acceptImproperIntersection=false) SIP_HOLDGIL
Compute the intersection between two segments.
double azimuth() const SIP_HOLDGIL
Returns the azimuth.
Definition: qgsellipse.h:139
static QgsCircle minimalCircleFrom3Points(const QgsPoint &pt1, const QgsPoint &pt2, const QgsPoint &pt3, double epsilon=1E-8) SIP_HOLDGIL
Constructs the smallest circle from 3 points.
Definition: qgscircle.cpp:223
A class to represent a 2D point.
Definition: qgspointxy.h:43
bool qgsDoubleNear(double a, double b, double epsilon=4 *std::numeric_limits< double >::epsilon())
Compare two doubles (but allow some difference)
Definition: qgis.h:316
int innerTangents(const QgsCircle &other, QgsPointXY &line1P1, QgsPointXY &line1P2, QgsPointXY &line2P1, QgsPointXY &line2P2) const
Calculates the inner tangent points between this circle and an other circle.
Definition: qgscircle.cpp:273
static QgsCircle fromCenterDiameter(const QgsPoint &center, double diameter, double azimuth=0) SIP_HOLDGIL
Constructs a circle by a center point and a diameter.
Definition: qgscircle.cpp:176
double semiMajorAxis() const SIP_HOLDGIL
Returns the semi-major axis.
Definition: qgsellipse.h:127
static QgsCircle from3Points(const QgsPoint &pt1, const QgsPoint &pt2, const QgsPoint &pt3, double epsilon=1E-8) SIP_HOLDGIL
Constructs a circle by 3 points on the circle.
Definition: qgscircle.cpp:83
static QgsPoint midpoint(const QgsPoint &pt1, const QgsPoint &pt2) SIP_HOLDGIL
Returns a middle point between points pt1 and pt2.
QVector< QgsPoint > northQuadrant() const
The four quadrants of the ellipse.
Definition: qgscircle.cpp:316
Triangle geometry type.
Definition: qgstriangle.h:33
static QgsCircle from2Points(const QgsPoint &pt1, const QgsPoint &pt2) SIP_HOLDGIL
Constructs a circle by 2 points on the circle.
Definition: qgscircle.cpp:37
bool tangentToPoint(const QgsPointXY &p, QgsPointXY &pt1, QgsPointXY &pt2) const
Calculates the tangent points between this circle and the point p.
Definition: qgscircle.cpp:262
QgsCircle inscribedCircle() const SIP_HOLDGIL
Inscribed circle of the triangle.
static double lineAngle(double x1, double y1, double x2, double y2) SIP_HOLDGIL
Calculates the direction of line joining two points in radians, clockwise from the north direction...
void setSemiMajorAxis(double semiMajorAxis) override SIP_HOLDGIL
Inherited method.
Definition: qgscircle.cpp:304
QgsPoint mCenter
Definition: qgsellipse.h:252
void setX(double x) SIP_HOLDGIL
Sets the point&#39;s x-coordinate.
Definition: qgspoint.h:269
double mSemiMajorAxis
Definition: qgsellipse.h:253
static bool setZValueFromPoints(const QgsPointSequence &points, QgsPoint &point)
A Z dimension is added to point if one of the point in the list points is in 3D.
double mAzimuth
Definition: qgsellipse.h:255
QString qgsDoubleToString(double a, int precision=17)
Returns a string representation of a double.
Definition: qgis.h:276
QString asWkt(int precision=17) const override
Returns a WKT representation of the geometry.
Definition: qgspoint.cpp:259
double perimeter() const override SIP_HOLDGIL
The circumference of the ellipse using first approximation of Ramanujan.
Definition: qgscircle.cpp:299
static Type dropZ(Type type) SIP_HOLDGIL
Drops the z dimension (if present) for a WKB type and returns the new type.
Definition: qgswkbtypes.h:1207
bool contains(const QgsPoint &point, double epsilon=1E-8) const
Returns true if the circle contains the point.
Definition: qgscircle.cpp:350
Point geometry type, with support for z-dimension and m-values.
Definition: qgspoint.h:37
double mSemiMinorAxis
Definition: qgsellipse.h:254
QVector< QgsPoint > QgsPointSequence
double area() const override SIP_HOLDGIL
The area of the ellipse.
Definition: qgscircle.cpp:294
static QgsCircle fromCenterPoint(const QgsPoint &center, const QgsPoint &pt1) SIP_HOLDGIL
Constructs a circle by a center point and another point.
Definition: qgscircle.cpp:181
static bool tangentPointAndCircle(const QgsPointXY &center, double radius, const QgsPointXY &p, QgsPointXY &pt1, QgsPointXY &pt2) SIP_HOLDGIL
Calculates the tangent points between the circle with the specified center and radius and the point p...
void setSemiMinorAxis(double semiMinorAxis) override SIP_HOLDGIL
Inherited method.
Definition: qgscircle.cpp:310
double distance(double x, double y) const SIP_HOLDGIL
Returns the Cartesian 2D distance between this point and a specified x, y coordinate.
Definition: qgspoint.h:332
virtual bool isEmpty() const SIP_HOLDGIL
An ellipse is empty if axes are equal to 0.
Definition: qgsellipse.cpp:118
bool is3D() const SIP_HOLDGIL
Returns true if the geometry is 3D and contains a z-value.
bool convertTo(QgsWkbTypes::Type type) override
Converts the geometry to a specified type.
Definition: qgspoint.cpp:610
QgsWkbTypes::Type wkbType() const SIP_HOLDGIL
Returns the WKB type of the geometry.
static int circleCircleIntersections(QgsPointXY center1, double radius1, QgsPointXY center2, double radius2, QgsPointXY &intersection1, QgsPointXY &intersection2) SIP_HOLDGIL
Calculates the intersections points between the circle with center center1 and radius radius1 and the...
double z
Definition: qgspoint.h:43
Circular string geometry type.
QgsRectangle boundingBox() const override
Returns the minimal bounding box for the ellipse.
Definition: qgscircle.cpp:355
QgsCircularString * toCircularString(bool oriented=false) const
Returns a circular string from the circle.
Definition: qgscircle.cpp:327
QgsPoint center() const SIP_HOLDGIL
Returns the center point.
Definition: qgsellipse.h:121
Ellipse geometry type.
Definition: qgsellipse.h:39
virtual QVector< QgsPoint > quadrant() const
The four quadrants of the ellipse.
Definition: qgsellipse.cpp:177
int intersections(const QgsCircle &other, QgsPoint &intersection1, QgsPoint &intersection2, bool useZ=false) const
Calculates the intersections points between this circle and an other circle.
Definition: qgscircle.cpp:239
int outerTangents(const QgsCircle &other, QgsPointXY &line1P1, QgsPointXY &line1P2, QgsPointXY &line2P1, QgsPointXY &line2P2) const
Calculates the outer tangent points between this circle and an other circle.
Definition: qgscircle.cpp:267
QString toString(int pointPrecision=17, int radiusPrecision=17, int azimuthPrecision=2) const override
returns a string representation of the ellipse.
Definition: qgscircle.cpp:360
double x
Definition: qgspoint.h:41
virtual QgsPointSequence points(unsigned int segments=36) const
Returns a list of points with segmentation from segments.
Definition: qgsellipse.cpp:188