30 const double t = ( ( ptX - x1 ) * dx + ( ptY - y1 ) * dy ) / ( dx * dx + dy * dy );
46 const double dist = dx * dx + dy * dy;
61 const double f1 = x - x1;
62 const double f2 = y2 - y1;
63 const double f3 = y - y1;
64 const double f4 = x2 - x1;
65 const double test = ( f1 * f2 - f3 * f4 );
67 return qgsDoubleNear( test, 0.0 ) ? 0 : ( test < 0 ? -1 : 1 );
70void QgsGeometryUtilsBase::pointOnLineWithDistance(
double x1,
double y1,
double x2,
double y2,
double distance,
double &x,
double &y,
double *z1,
double *z2,
double *z,
double *m1,
double *m2,
double *m )
72 const double dx = x2 - x1;
73 const double dy = y2 - y1;
74 const double length = std::sqrt( dx * dx + dy * dy );
87 const double scaleFactor = distance / length;
88 x = x1 + dx * scaleFactor;
89 y = y1 + dy * scaleFactor;
91 *z = *z1 + ( *z2 - *z1 ) * scaleFactor;
93 *m = *m1 + ( *m2 - *m1 ) * scaleFactor;
100 const double mX = x1 + ( x2 - x1 ) * proportion;
101 const double mY = y1 + ( y2 - y1 ) * proportion;
102 const double pX = x1 - x2;
103 const double pY = y1 - y2;
104 double normalX = -pY;
106 const double normalLength = sqrt( ( normalX * normalX ) + ( normalY * normalY ) );
107 normalX /= normalLength;
108 normalY /= normalLength;
110 *x = mX + offset * normalX;
111 *y = mY + offset * normalY;
116 const double angle = std::atan2( dy, dx ) * 180 / M_PI;
121 else if ( angle > 360 )
130 if ( angle3 >= angle1 )
132 return !( angle2 > angle1 && angle2 < angle3 );
136 return !( angle2 > angle1 || angle2 < angle3 );
144 if ( angle2 < angle1 )
146 return ( angle <= angle1 && angle >= angle2 );
150 return ( angle <= angle1 || angle >= angle2 );
155 if ( angle2 > angle1 )
157 return ( angle >= angle1 && angle <= angle2 );
161 return ( angle >= angle1 || angle <= angle2 );
174 double dx21, dy21, dx31, dy31, h21, h31, d;
179 centerX = ( x1 + x2 ) / 2.0;
180 centerY = ( y1 + y2 ) / 2.0;
181 radius = std::sqrt( std::pow( centerX - x1, 2.0 ) + std::pow( centerY - y1, 2.0 ) );
191 h21 = std::pow( dx21, 2.0 ) + std::pow( dy21, 2.0 );
192 h31 = std::pow( dx31, 2.0 ) + std::pow( dy31, 2.0 );
195 d = 2 * ( dx21 * dy31 - dx31 * dy21 );
205 centerX = x1 + ( h21 * dy31 - h31 * dy21 ) / d;
206 centerY = y1 - ( h21 * dx31 - h31 * dx21 ) / d;
207 radius = std::sqrt( std::pow( centerX - x1, 2.0 ) + std::pow( centerY - y1, 2.0 ) );
216 double length = M_PI / 180.0 * radius *
sweepAngle( centerX, centerY, x1, y1, x2, y2, x3, y3 );
230 if ( p3Angle >= p1Angle )
232 if ( p2Angle > p1Angle && p2Angle < p3Angle )
234 return ( p3Angle - p1Angle );
238 return ( - ( p1Angle + ( 360 - p3Angle ) ) );
243 if ( p2Angle < p1Angle && p2Angle > p3Angle )
245 return ( -( p1Angle - p3Angle ) );
249 return ( p3Angle + ( 360 - p1Angle ) );
260 return zm1 + ( zm2 - zm1 ) * ( angle - a1 ) / ( a2 - a1 );
262 return zm2 + ( zm3 - zm2 ) * ( angle - a2 ) / ( a3 - a2 );
268 return zm1 + ( zm2 - zm1 ) * ( a1 - angle ) / ( a1 - a2 );
270 return zm2 + ( zm3 - zm2 ) * ( a2 - angle ) / ( a2 - a3 );
275 double clippedAngle = angle;
276 if ( clippedAngle >= M_PI * 2 || clippedAngle <= -2 * M_PI )
278 clippedAngle = std::fmod( clippedAngle, 2 * M_PI );
280 if ( clippedAngle < 0.0 )
282 clippedAngle += 2 * M_PI;
291 if ( x <= left && y <= bottom )
293 const double dx = left - x;
294 const double dy = bottom - y;
302 else if ( x >= right && y >= top )
304 const double dx = x - right;
305 const double dy = y - top;
313 else if ( x >= right && y <= bottom )
315 const double dx = x - right;
316 const double dy = bottom - y;
324 else if ( x <= left && y >= top )
326 const double dx = left - x;
327 const double dy = y - top;
335 else if ( x <= left )
337 else if ( x >= right )
339 else if ( y <= bottom )
345 const double smallestX = std::min( right - x, x - left );
346 const double smallestY = std::min( top - y, y - bottom );
347 if ( smallestX < smallestY )
350 if ( right - x < x - left )
358 if ( top - y < y - bottom )
365void QgsGeometryUtilsBase::perpendicularCenterSegment(
double pointx,
double pointy,
double segmentPoint1x,
double segmentPoint1y,
double segmentPoint2x,
double segmentPoint2y,
double &perpendicularSegmentPoint1x,
double &perpendicularSegmentPoint1y,
double &perpendicularSegmentPoint2x,
double &perpendicularSegmentPoint2y,
double desiredSegmentLength )
367 QgsVector segmentVector =
QgsVector( segmentPoint2x - segmentPoint1x, segmentPoint2y - segmentPoint1y );
369 if ( desiredSegmentLength != 0 )
371 perpendicularVector = perpendicularVector.
normalized() * ( desiredSegmentLength ) / 2;
374 perpendicularSegmentPoint1x = pointx - perpendicularVector.
x();
375 perpendicularSegmentPoint1y = pointy - perpendicularVector.
y();
376 perpendicularSegmentPoint2x = pointx + perpendicularVector.
x();
377 perpendicularSegmentPoint2y = pointy + perpendicularVector.
y();
382 const double at = std::atan2( y2 - y1, x2 - x1 );
383 const double a = -at + M_PI_2;
389 const double angle1 = std::atan2( y1 - y2, x1 - x2 );
390 const double angle2 = std::atan2( y3 - y2, x3 - x2 );
404 const double a1 =
lineAngle( x1, y1, x2, y2 );
405 const double a2 =
lineAngle( x2, y2, x3, y3 );
413 double clockwiseDiff = 0.0;
416 clockwiseDiff = a2 - a1;
420 clockwiseDiff = a2 + ( 2 * M_PI - a1 );
422 const double counterClockwiseDiff = 2 * M_PI - clockwiseDiff;
424 double resultAngle = 0;
425 if ( clockwiseDiff <= counterClockwiseDiff )
427 resultAngle = a1 + clockwiseDiff / 2.0;
431 resultAngle = a1 - counterClockwiseDiff / 2.0;
442 if ( u3.
length() == 0 )
return 1;
459 if ( std::fabs( u3.
x() ) <= epsilon &&
460 std::fabs( u3.
y() ) <= epsilon &&
461 std::fabs( u3.
z() ) <= epsilon )
473 if ( !( std::fabs( b1 ) > epsilon ) )
478 if ( !( a2 != -1 && a2 != 1 ) )
484 const double r1 = ( c2 - b2 * c1 / b1 ) / ( a2 - b2 * a1 / b1 );
492 const double d = v1.
y() * v2.
x() - v1.
x() * v2.
y();
497 const double dx = p2x - p1x;
498 const double dy = p2y - p1y;
499 const double k = ( dy * v2.
x() - dx * v2.
y() ) / d;
501 intersectionX = p1x + v1.
x() * k;
502 intersectionY = p1y + v1.
y() * k;
507static bool equals(
double p1x,
double p1y,
double p2x,
double p2y,
double epsilon = 1e-8 )
512bool QgsGeometryUtilsBase::segmentIntersection(
double p1x,
double p1y,
double p2x,
double p2y,
double q1x,
double q1y,
double q2x,
double q2y,
double &intersectionPointX,
double &intersectionPointY,
bool &isIntersection,
double tolerance,
bool acceptImproperIntersection )
514 isIntersection =
false;
515 intersectionPointX = intersectionPointY = std::numeric_limits<double>::quiet_NaN();
519 const double vl = v.
length();
520 const double wl = w.
length();
529 if ( !
lineIntersection( p1x, p1y, v, q1x, q1y, w, intersectionPointX, intersectionPointY ) )
534 isIntersection =
true;
535 if ( acceptImproperIntersection )
537 if ( ( equals( p1x, p1y, q1x, q1y ) ) || ( equals( p1x, p1y, q2x, q2y ) ) )
539 intersectionPointX = p1x;
540 intersectionPointY = p1y;
543 else if ( ( equals( p1x, p1y, q2x, q2y ) ) || ( equals( p2x, p2y, q2x, q2y ) ) )
545 intersectionPointX = p2x;
546 intersectionPointY = p2y;
553 qgsDoubleNear(
sqrDistToLine( p1x, p1y, q1x, q1y, q2x, q2y, x, y, tolerance ), 0.0, tolerance ) ||
555 qgsDoubleNear(
sqrDistToLine( p2x, p2y, q1x, q1y, q2x, q2y, x, y, tolerance ), 0.0, tolerance ) ||
557 qgsDoubleNear(
sqrDistToLine( q1x, q1y, p1x, p1y, p2x, p2y, x, y, tolerance ), 0.0, tolerance ) ||
559 qgsDoubleNear(
sqrDistToLine( q2x, q2y, p1x, p1y, p2x, p2y, x, y, tolerance ), 0.0, tolerance )
566 const double lambdav =
QgsVector( intersectionPointX - p1x, intersectionPointY - p1y ) * v;
567 if ( lambdav < 0. + tolerance || lambdav > vl - tolerance )
570 const double lambdaw =
QgsVector( intersectionPointX - q1x, intersectionPointY - q1y ) * w;
571 return !( lambdaw < 0. + tolerance || lambdaw >= wl - tolerance );
585 double ptInterX = 0.0, ptInterY = 0.0;
586 bool isIntersection =
false;
595 intersection.
set( ptInterX, ptInterY, La1.
z() );
635 if ( !firstIsDone || !secondIsDone )
641 intersection = ( X1 + X2 ) / 2.0;
647 return 0.5 * std::abs( ( aX - cX ) * ( bY - aY ) - ( aX - bX ) * ( cY - aY ) );
652 const double dxp = px - x1;
653 const double dyp = py - y1;
655 const double dxl = x2 - x1;
656 const double dyl = y2 - y1;
658 return std::sqrt( ( dxp * dxp ) + ( dyp * dyp ) ) / std::sqrt( ( dxl * dxl ) + ( dyl * dyl ) );
662 double weightB,
double weightC,
double &pointX,
double &pointY )
665 if ( weightB + weightC > 1 )
667 weightB = 1 - weightB;
668 weightC = 1 - weightC;
671 const double rBx = weightB * ( bX - aX );
672 const double rBy = weightB * ( bY - aY );
673 const double rCx = weightC * ( cX - aX );
674 const double rCy = weightC * ( cY - aY );
676 pointX = rBx + rCx + aX;
677 pointY = rBy + rCy + aY;
682 return qgsDoubleNear( x1 * ( y2 - y3 ) + x2 * ( y3 - y1 ) + x3 * ( y1 - y2 ), 0, epsilon );
687 const double dx = x2 - x1;
688 const double dy = y2 - y1;
689 return ( std::atan2( dx, dy ) * 180.0 / M_PI );
693 double &pointX,
double &pointY,
double &angle )
695 angle = (
azimuth( aX, aY, bX, bY ) +
azimuth( cX, cY, dX, dY ) ) / 2.0;
697 bool intersection =
false;
698 QgsGeometryUtilsBase::segmentIntersection( aX, aY, bX, bY, cX, cY, dX, dY, pointX, pointY, intersection );
703void QgsGeometryUtilsBase::project(
double aX,
double aY,
double aZ,
double distance,
double azimuth,
double inclination,
double &resultX,
double &resultY,
double &resultZ )
705 const double radsXy =
azimuth * M_PI / 180.0;
706 double dx = 0.0, dy = 0.0, dz = 0.0;
708 inclination = std::fmod( inclination, 360.0 );
710 if ( std::isnan( aZ ) &&
qgsDoubleNear( inclination, 90.0 ) )
712 dx = distance * std::sin( radsXy );
713 dy = distance * std::cos( radsXy );
717 const double radsZ = inclination * M_PI / 180.0;
718 dx = distance * std::sin( radsZ ) * std::sin( radsXy );
719 dy = distance * std::sin( radsZ ) * std::cos( radsXy );
720 dz = distance * std::cos( radsZ );
729 double &pointX,
double &pointY )
731 const double angle = (
azimuth( aX, aY, bX, bY ) +
azimuth( aX, aY, cX, cY ) ) / 2.0;
733 bool intersection =
false;
734 double dX = 0.0, dY = 0.0, dZ = 0.0;
735 project( aX, aY, std::numeric_limits<double>::quiet_NaN(), 1.0, angle, 90.0, dX, dY, dZ );
736 segmentIntersection( bX, bY, cX, cY, aX, aY, dX, dY, pointX, pointY, intersection );
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.
static void pointOnLineWithDistance(double x1, double y1, double x2, double y2, double distance, double &x, double &y, double *z1=nullptr, double *z2=nullptr, double *z=nullptr, double *m1=nullptr, double *m2=nullptr, double *m=nullptr)
Calculates the point a specified distance from (x1, y1) toward a second point (x2,...
static double angleBetweenThreePoints(double x1, double y1, double x2, double y2, double x3, double y3)
Calculates the angle between the lines AB and BC, where AB and BC described by points a,...
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...
static bool circleClockwise(double angle1, double angle2, double angle3)
Returns true if the circle defined by three angles is ordered clockwise.
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.
static double linePerpendicularAngle(double x1, double y1, double x2, double y2)
Calculates the perpendicular angle to a line joining two points.
static bool bisector(double aX, double aY, double bX, double bY, double cX, double cY, double &pointX, double &pointY)
Returns the point (pointX, pointY) forming the bisector from point (aX, aY) to the segment (bX,...
static double lineAngle(double x1, double y1, double x2, double y2)
Calculates the direction of line joining two points in radians, clockwise from the north direction.
static double averageAngle(double x1, double y1, double x2, double y2, double x3, double y3)
Calculates the average angle (in radians) between the two linear segments from (x1,...
static bool skewLinesProjection(const QgsVector3D &P1, const QgsVector3D &P12, const QgsVector3D &P2, const QgsVector3D &P22, QgsVector3D &X1, double epsilon=0.0001)
A method to project one skew line onto another.
static int closestSideOfRectangle(double right, double bottom, double left, double top, double x, double y)
Returns a number representing the closest side of a rectangle defined by /a right,...
static void project(double aX, double aY, double aZ, double distance, double azimuth, double inclination, double &resultX, double &resultY, double &resultZ)
Returns coordinates of a point which corresponds to this point projected by a specified distance with...
static double interpolateArcValue(double angle, double a1, double a2, double a3, double zm1, double zm2, double zm3)
Interpolate a value at given angle on circular arc given values (zm1, zm2, zm3) at three different an...
static void perpendicularOffsetPointAlongSegment(double x1, double y1, double x2, double y2, double proportion, double offset, double *x, double *y)
Calculates a point a certain proportion of the way along the segment from (x1, y1) to (x2,...
static double skewLinesDistance(const QgsVector3D &P1, const QgsVector3D &P12, const QgsVector3D &P2, const QgsVector3D &P22)
An algorithm to calculate the shortest distance between two skew lines.
static double pointFractionAlongLine(double x1, double y1, double x2, double y2, double px, double py)
Given the line (x1, y1) to (x2, y2) and a point (px, py) returns the fraction of the line length at w...
static double normalizedAngle(double angle)
Ensures that an angle is in the range 0 <= angle < 2 pi.
static void weightedPointInTriangle(double aX, double aY, double bX, double bY, double cX, double cY, double weightB, double weightC, double &pointX, double &pointY)
Returns a weighted point inside the triangle denoted by the points (aX, aY), (bX, bY) and (cX,...
static void perpendicularCenterSegment(double centerPointX, double centerPointY, double segmentPoint1x, double segmentPoint1y, double segmentPoint2x, double segmentPoint2y, double &perpendicularSegmentPoint1x, double &perpendicularSegmentPoint1y, double &perpendicularSegmentPoint2x, double &perpendicularSegmentPoint2y, double segmentLength=0)
Create a perpendicular line segment to a given segment [segmentPoint1,segmentPoint2] with its center ...
static bool circleAngleBetween(double angle, double angle1, double angle2, bool clockwise)
Returns true if, in a circle, angle is between angle1 and angle2.
static double triangleArea(double aX, double aY, double bX, double bY, double cX, double cY)
Returns the area of the triangle denoted by the points (aX, aY), (bX, bY) and (cX,...
static bool angleBisector(double aX, double aY, double bX, double bY, double cX, double cY, double dX, double dY, double &pointX, double &pointY, double &angle)
Returns the point (pointX, pointY) forming the bisector from segment (aX aY) (bX bY) and segment (bX,...
static double sqrDistToLine(double ptX, double ptY, double x1, double y1, double x2, double y2, double &minDistX, double &minDistY, double epsilon)
Returns the squared distance between a point and a line.
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,...
static bool segmentIntersection(double p1x, double p1y, double p2x, double p2y, double q1x, double q1y, double q2x, double q2y, double &intersectionPointX, double &intersectionPointY, bool &isIntersection, double tolerance=1e-8, bool acceptImproperIntersection=false)
Compute the intersection between two segments.
static int leftOfLine(const double x, const double y, const double x1, const double y1, const double x2, const double y2)
Returns a value < 0 if the point (x, y) is left of the line from (x1, y1) -> (x2, y2).
static bool linesIntersection3D(const QgsVector3D &La1, const QgsVector3D &La2, const QgsVector3D &Lb1, const QgsVector3D &Lb2, QgsVector3D &intersection)
An algorithm to calculate an (approximate) intersection of two lines in 3D.
static bool pointsAreCollinear(double x1, double y1, double x2, double y2, double x3, double y3, double epsilon)
Given the points (x1, y1), (x2, y2) and (x3, y3) returns true if these points can be considered colli...
static void circleCenterRadius(double x1, double y1, double x2, double y2, double x3, double y3, double &radius, double ¢erX, double ¢erY)
Returns radius and center of the circle through (x1 y1), (x2 y2), (x3 y3)
static bool lineIntersection(double p1x, double p1y, QgsVector v1, double p2x, double p2y, QgsVector v2, double &intersectionX, double &intersectionY)
Computes the intersection between two lines.
static double azimuth(double x1, double y1, double x2, double y2)
Calculates Cartesian azimuth between points (x1, y1) and (x2, y2) (clockwise in degree,...
Class for storage of 3D vectors similar to QVector3D, with the difference that it uses double precisi...
double y() const
Returns Y coordinate.
double z() const
Returns Z coordinate.
static double dotProduct(const QgsVector3D &v1, const QgsVector3D &v2)
Returns the dot product of two vectors.
double x() const
Returns X coordinate.
void normalize()
Normalizes the current vector in place.
static QgsVector3D crossProduct(const QgsVector3D &v1, const QgsVector3D &v2)
Returns the cross product of two vectors.
void set(double x, double y, double z)
Sets vector coordinates.
double length() const
Returns the length of the vector.
A class to represent a vector.
double y() const
Returns the vector's y-component.
QgsVector normalized() const
Returns the vector's normalized (or "unit") vector (ie same angle but length of 1....
QgsVector perpVector() const
Returns the perpendicular vector to this vector (rotated 90 degrees counter-clockwise)
double x() const
Returns the vector's x-component.
double length() const
Returns the length of the vector.
bool qgsDoubleNear(double a, double b, double epsilon=4 *std::numeric_limits< double >::epsilon())
Compare two doubles (but allow some difference)