Skip to content
Permalink
Browse files

Merge pull request #5858 from lbartoletti/segment_intersection

Segment intersection
  • Loading branch information
nyalldawson committed Dec 22, 2017
2 parents 36c4ac4 + 313417d commit 4ff72de62e3e95ab8bbd63d68964c0f2809311f5
@@ -1465,7 +1465,8 @@ QgsGeometryUtils {#qgis_api_break_3_0_QgsGeometryUtils}

- componentType enum has been renamed to ComponentType and its members were CamelCased too: VERTEX, RING and PART become Vertex, Ring and Part, respectively.
- adjacentVertices was removed - use QgsAbstractGeometry.adjacentVertices instead.

- segmentIntersection takes an optional boolean parameter "acceptImproperIntersection" returning true even if the intersection is improper.
Takes another boolean argument "isIntersection" returning if there is an intersection or not. The "inter" parameter has been renamed "intersectionPoint".

QgsGPSConnectionRegistry {#qgis_api_break_3_0_QgsGPSConnectionRegistry}
------------------------
@@ -100,18 +100,42 @@ Returns the squared distance between a point and a line.
:return: Whether the lines intersect
%End

static bool segmentIntersection( const QgsPoint &p1, const QgsPoint &p2, const QgsPoint &q1, const QgsPoint &q2, QgsPoint &inter /Out/, double tolerance );
static bool segmentIntersection( const QgsPoint &p1, const QgsPoint &p2, const QgsPoint &q1, const QgsPoint &q2, QgsPoint &intersectionPoint /Out/, bool &isIntersection /Out/, const double tolerance = 1e-8, bool acceptImproperIntersection = false );
%Docstring
Compute the intersection between two segments

:param p1: First segment start point
:param p2: First segment end point
:param q1: Second segment start point
:param q2: Second segment end point
:param inter: Output parameter, the intersection point
:param intersectionPoint: Output parameter, the intersection point
:param isIntersection: Output parameter, return true if an intersection is found
:param tolerance: The tolerance to use
:param acceptImproperIntersection: By default, this method returns true only if segments have proper intersection. If set true, returns also true if segments have improper intersection (end of one segment on other segment ; continuous segments).

:return: Whether the segments intersect
* Example:
\code{.py}
ret = QgsGeometryUtils.segmentIntersection( QgsPoint( 0, 0 ), QgsPoint( 0, 1 ), QgsPoint( 1, 1 ), QgsPoint( 1, 0 ) )
ret[0], ret[1].asWkt(), ret[2]
# Whether the segments intersect, the intersection point, is intersect
# (False, 'Point (0 0)', False)
ret = QgsGeometryUtils.segmentIntersection( QgsPoint( 0, 0 ), QgsPoint( 0, 5 ), QgsPoint( 0, 5 ), QgsPoint( 1, 5 ) )
ret[0], ret[1].asWkt(), ret[2]
# (False, 'Point (0 5)', True)
ret = QgsGeometryUtils.segmentIntersection( QgsPoint( 0, 0 ), QgsPoint( 0, 5 ), QgsPoint( 0, 5 ), QgsPoint( 1, 5 ), acceptImproperIntersection=True )
ret[0], ret[1].asWkt(), ret[2]
# (True, 'Point (0 5)', True)
ret = QgsGeometryUtils.segmentIntersection( QgsPoint( 0, 0 ), QgsPoint( 0, 5 ), QgsPoint( 0, 2 ), QgsPoint( 1, 5 ) )
ret[0], ret[1].asWkt(), ret[2]
# (False, 'Point (0 2)', True)
ret = QgsGeometryUtils.segmentIntersection( QgsPoint( 0, 0 ), QgsPoint( 0, 5 ), QgsPoint( 0, 2 ), QgsPoint( 1, 5 ), acceptImproperIntersection=True )
ret[0], ret[1].asWkt(), ret[2]
# (True, 'Point (0 2)', True)
ret = QgsGeometryUtils.segmentIntersection( QgsPoint( 0, -5 ), QgsPoint( 0, 5 ), QgsPoint( 2, 0 ), QgsPoint( -1, 0 ) )
ret[0], ret[1].asWkt(), ret[2]
# (True, 'Point (0 0)', True)
\endcode
%End

static QgsPoint projPointOnSegment( const QgsPoint &p, const QgsPoint &s1, const QgsPoint &s2 );
@@ -257,6 +257,7 @@ namespace QgsGeometryCheckerUtils
{
QList<QgsPoint> intersections;
QgsPoint inter;
bool intersection = false;
for ( int i = 0, n = line1->vertexCount() - 1; i < n; ++i )
{
for ( int j = 0, m = line2->vertexCount() - 1; j < m; ++j )
@@ -265,7 +266,7 @@ namespace QgsGeometryCheckerUtils
QgsPoint p2 = line1->vertexAt( QgsVertexId( 0, 0, i + 1 ) );
QgsPoint q1 = line2->vertexAt( QgsVertexId( 0, 0, j ) );
QgsPoint q2 = line2->vertexAt( QgsVertexId( 0, 0, j + 1 ) );
if ( QgsGeometryUtils::segmentIntersection( p1, p2, q1, q2, inter, tol ) )
if ( QgsGeometryUtils::segmentIntersection( p1, p2, q1, q2, inter, intersection, tol ) )
{
intersections.append( inter );
}
@@ -114,7 +114,8 @@ void QgsGeometrySelfIntersectionCheck::fixError( QgsGeometryCheckError *error, i
QgsPoint p2 = geom->vertexAt( QgsVertexId( vidx.part, vidx.ring, ( inter.segment1 + 1 ) % nVerts ) );
QgsPoint q2 = geom->vertexAt( QgsVertexId( vidx.part, vidx.ring, ( inter.segment2 + 1 ) % nVerts ) );
QgsPoint s;
if ( !QgsGeometryUtils::segmentIntersection( p1, p2, q1, q2, s, mContext->tolerance ) )
bool intersection = false;
if ( !QgsGeometryUtils::segmentIntersection( p1, p2, q1, q2, s, intersection, mContext->tolerance ) )
{
error->setObsolete();
return;
@@ -25,6 +25,7 @@
#include "qgslinestring.h"
#include "qgsmultipolygon.h"
#include "qgsspinbox.h"
#include "qgsgeometryutils.h"
#include <memory>
#include <QMouseEvent>

@@ -61,9 +62,11 @@ void QgsMapToolCircle2TangentsPoint::cadCanvasReleaseEvent( QgsMapMouseEvent *e
}
if ( mPoints.size() == 4 )
{
QgsPointXY ptInter = intersect( QgsPointXY( mPoints.at( 0 ) ), QgsPointXY( mPoints.at( 1 ) ),
QgsPointXY( mPoints.at( 2 ) ), QgsPointXY( mPoints.at( 3 ) ) );
if ( ptInter == QgsPointXY() )
bool isIntersect = false;
QgsPoint ptInter;
QgsGeometryUtils::segmentIntersection( mPoints.at( 0 ), mPoints.at( 1 ),
mPoints.at( 2 ), mPoints.at( 3 ), ptInter, isIntersect );
if ( !isIntersect )
{
QgisApp::instance()->messageBar()->pushMessage( tr( "Error" ), tr( "Segments are parallels" ),
QgsMessageBar::CRITICAL, QgisApp::instance()->messageTimeout() );
@@ -144,64 +147,6 @@ void QgsMapToolCircle2TangentsPoint::cadCanvasMoveEvent( QgsMapMouseEvent *e )
}
}

QgsPointXY QgsMapToolCircle2TangentsPoint::intersect( QgsPointXY seg1_pt1, QgsPointXY seg1_pt2, QgsPointXY seg2_pt1, QgsPointXY seg2_pt2 )
{
/*
* Public domain function by Darel Rex Finley, 2006
* http://alienryderflex.com/intersect/
*/
QgsPointXY ptInter;

double Ax = seg1_pt1.x();
double Ay = seg1_pt1.y();
double Bx = seg1_pt2.x();
double By = seg1_pt2.y();

double Cx = seg2_pt1.x();
double Cy = seg2_pt1.y();
double Dx = seg2_pt2.x();
double Dy = seg2_pt2.y();

if ( ( ( Ax == Bx ) && ( Ay == By ) ) || ( ( Cx == Dx ) && ( Cy == Dy ) ) )
return ptInter;

// (1) Translate the system so that point A is on the origin.
Bx -= Ax;
By -= Ay;
Cx -= Ax;
Cy -= Ay;
Dx -= Ax;
Dy -= Ay;

// Discover the length of segment A-B
double distAB = sqrt( Bx * Bx + By * By );

// (2) Rotate the system so that point B is on the positive X axis.
double theCos = Bx / distAB;
double theSin = By / distAB;
double newX = Cx * theCos + Cy * theSin;
Cy = Cy * theCos - Cx * theSin;
Cx = newX;
newX = Dx * theCos + Dy * theSin;
Dy = Dy * theCos - Dx * theSin;
Dx = newX;

// Fail if the lines are parallel.
if ( Cy == Dy )
return ptInter;

// (3) Discover the position of the intersection point along line A-B.
double ABpos = Dx + ( Cx - Dx ) * Dy / ( Dy - Cy );

// (4) Apply the discovered position to line A-B
// in the original coordinate system.
ptInter.setX( Ax + ABpos * theCos );
ptInter.setY( Ay + ABpos * theSin );

// Success
return ptInter;
}

void QgsMapToolCircle2TangentsPoint::getPossibleCenter( )
{

@@ -226,20 +171,20 @@ void QgsMapToolCircle2TangentsPoint::getPossibleCenter( )
QgsGeometry line2m = line2.offsetCurve( - mRadius, 8, QgsGeometry::JoinStyleBevel, 5 );
QgsGeometry line2p = line2.offsetCurve( + mRadius, 8, QgsGeometry::JoinStyleBevel, 5 );

QgsPointXY p1 = intersect( line1m.asPolyline().at( 0 ), line1m.asPolyline().at( 1 ),
line2m.asPolyline().at( 0 ), line2m.asPolyline().at( 1 ) );
QgsPointXY p2 = intersect( line1m.asPolyline().at( 0 ), line1m.asPolyline().at( 1 ),
line2p.asPolyline().at( 0 ), line2p.asPolyline().at( 1 ) );
QgsPointXY p3 = intersect( line1p.asPolyline().at( 0 ), line1p.asPolyline().at( 1 ),
line2m.asPolyline().at( 0 ), line2m.asPolyline().at( 1 ) );
QgsPointXY p4 = intersect( line1p.asPolyline().at( 0 ), line1p.asPolyline().at( 1 ),
line2p.asPolyline().at( 0 ), line2p.asPolyline().at( 1 ) );

mCenters.append( p1 );
mCenters.append( p2 );
mCenters.append( p3 );
mCenters.append( p4 );

bool isIntersect = false;
QgsPoint inter;
QgsGeometryUtils::segmentIntersection( QgsPoint( line1m.asPolyline().at( 0 ) ), QgsPoint( line1m.asPolyline().at( 1 ) ),
QgsPoint( line2m.asPolyline().at( 0 ) ), QgsPoint( line2m.asPolyline().at( 1 ) ), inter, isIntersect );
mCenters.append( QgsPointXY( inter ) );
QgsGeometryUtils::segmentIntersection( QgsPoint( line1m.asPolyline().at( 0 ) ), QgsPoint( line1m.asPolyline().at( 1 ) ),
QgsPoint( line2p.asPolyline().at( 0 ) ), QgsPoint( line2p.asPolyline().at( 1 ) ), inter, isIntersect );
mCenters.append( QgsPointXY( inter ) );
QgsGeometryUtils::segmentIntersection( QgsPoint( line1p.asPolyline().at( 0 ) ), QgsPoint( line1p.asPolyline().at( 1 ) ),
QgsPoint( line2m.asPolyline().at( 0 ) ), QgsPoint( line2m.asPolyline().at( 1 ) ), inter, isIntersect );
mCenters.append( QgsPointXY( inter ) );
QgsGeometryUtils::segmentIntersection( QgsPoint( line1p.asPolyline().at( 0 ) ), QgsPoint( line1p.asPolyline().at( 1 ) ),
QgsPoint( line2p.asPolyline().at( 0 ) ), QgsPoint( line2p.asPolyline().at( 1 ) ), inter, isIntersect );
mCenters.append( QgsPointXY( inter ) );
}
}

@@ -38,9 +38,6 @@ class QgsMapToolCircle2TangentsPoint: public QgsMapToolAddCircle
void radiusSpinBoxChanged( int radius );

private:
//! Return the point where segments are intersected. Method from QgsGeometryUtils doesn't work for special cases used by this tool.
QgsPointXY intersect( QgsPointXY seg1_pt1, QgsPointXY seg1_pt2, QgsPointXY seg2_pt1, QgsPointXY seg2_pt2 );

//! Compute 4 possible centers
void getPossibleCenter();

@@ -182,9 +182,16 @@ QgsCircle QgsCircle::fromCenterPoint( const QgsPoint &center, const QgsPoint &pt
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 )
{
QgsPoint p1, p2, p3;
QgsGeometryUtils::segmentIntersection( pt1_tg1, pt2_tg1, pt1_tg2, pt2_tg2, p1, epsilon );
QgsGeometryUtils::segmentIntersection( pt1_tg1, pt2_tg1, pt1_tg3, pt2_tg3, p2, epsilon );
QgsGeometryUtils::segmentIntersection( pt1_tg2, pt2_tg2, pt1_tg3, pt2_tg3, p3, epsilon );
bool isIntersect = false;
QgsGeometryUtils::segmentIntersection( pt1_tg1, pt2_tg1, pt1_tg2, pt2_tg2, p1, isIntersect, epsilon );
if ( !isIntersect )
return QgsCircle();
QgsGeometryUtils::segmentIntersection( pt1_tg1, pt2_tg1, pt1_tg3, pt2_tg3, p2, isIntersect, epsilon );
if ( !isIntersect )
return QgsCircle();
QgsGeometryUtils::segmentIntersection( pt1_tg2, pt2_tg2, pt1_tg3, pt2_tg3, p3, isIntersect, epsilon );
if ( !isIntersect )
return QgsCircle();

return QgsTriangle( p1, p2, p3 ).inscribedCircle();
}
@@ -251,28 +251,62 @@ bool QgsGeometryUtils::lineIntersection( const QgsPoint &p1, QgsVector v1, const
return true;
}

bool QgsGeometryUtils::segmentIntersection( const QgsPoint &p1, const QgsPoint &p2, const QgsPoint &q1, const QgsPoint &q2, QgsPoint &inter, double tolerance )
bool QgsGeometryUtils::segmentIntersection( const QgsPoint &p1, const QgsPoint &p2, const QgsPoint &q1, const QgsPoint &q2, QgsPoint &intersectionPoint, bool &isIntersection, const double tolerance, bool acceptImproperIntersection )
{
isIntersection = false;

QgsVector v( p2.x() - p1.x(), p2.y() - p1.y() );
QgsVector w( q2.x() - q1.x(), q2.y() - q1.y() );
double vl = v.length();
double wl = w.length();

if ( qgsDoubleNear( vl, 0, 0.000000000001 ) || qgsDoubleNear( wl, 0, 0.000000000001 ) )
if ( qgsDoubleNear( vl, 0.0, tolerance ) || qgsDoubleNear( wl, 0.0, tolerance ) )
{
return false;
}
v = v / vl;
w = w / wl;

if ( !QgsGeometryUtils::lineIntersection( p1, v, q1, w, inter ) )
if ( !QgsGeometryUtils::lineIntersection( p1, v, q1, w, intersectionPoint ) )
{
return false;
}

isIntersection = true;
if ( acceptImproperIntersection )
{
if ( ( p1 == q1 ) || ( p1 == q2 ) )
{
intersectionPoint = p1;
return true;
}
else if ( ( p2 == q1 ) || ( p2 == q2 ) )
{
intersectionPoint = p2;
return true;
}

double x, y;
if (
// intersectionPoint = p1
qgsDoubleNear( QgsGeometryUtils::sqrDistToLine( p1.x(), p1.y(), q1.x(), q1.y(), q2.x(), q2.y(), x, y, tolerance ), 0.0, tolerance ) ||
// intersectionPoint = p2
qgsDoubleNear( QgsGeometryUtils::sqrDistToLine( p2.x(), p2.y(), q1.x(), q1.y(), q2.x(), q2.y(), x, y, tolerance ), 0.0, tolerance ) ||
// intersectionPoint = q1
qgsDoubleNear( QgsGeometryUtils::sqrDistToLine( q1.x(), q1.y(), p1.x(), p1.y(), p2.x(), p2.y(), x, y, tolerance ), 0.0, tolerance ) ||
// intersectionPoint = q2
qgsDoubleNear( QgsGeometryUtils::sqrDistToLine( q2.x(), q2.y(), p1.x(), p1.y(), p2.x(), p2.y(), x, y, tolerance ), 0.0, tolerance )
)
{
return true;
}
}

double lambdav = QgsVector( inter.x() - p1.x(), inter.y() - p1.y() ) * v;
double lambdav = QgsVector( intersectionPoint.x() - p1.x(), intersectionPoint.y() - p1.y() ) * v;
if ( lambdav < 0. + tolerance || lambdav > vl - tolerance )
return false;

double lambdaw = QgsVector( inter.x() - q1.x(), inter.y() - q1.y() ) * w;
double lambdaw = QgsVector( intersectionPoint.x() - q1.x(), intersectionPoint.y() - q1.y() ) * w;
return !( lambdaw < 0. + tolerance || lambdaw >= wl - tolerance );
}

@@ -299,7 +333,8 @@ QVector<QgsGeometryUtils::SelfIntersection> QgsGeometryUtils::getSelfIntersectio
QgsPoint pl = geom->vertexAt( QgsVertexId( part, ring, l ) );

QgsPoint inter;
if ( !QgsGeometryUtils::segmentIntersection( pi, pj, pk, pl, inter, tolerance ) ) continue;
bool intersection = false;
if ( !QgsGeometryUtils::segmentIntersection( pi, pj, pk, pl, inter, intersection, tolerance ) ) continue;

SelfIntersection s;
s.segment1 = i;
@@ -106,11 +106,35 @@ class CORE_EXPORT QgsGeometryUtils
* \param p2 First segment end point
* \param q1 Second segment start point
* \param q2 Second segment end point
* \param inter Output parameter, the intersection point
* \param intersectionPoint Output parameter, the intersection point
* \param isIntersection Output parameter, return true if an intersection is found
* \param tolerance The tolerance to use
* \param acceptImproperIntersection By default, this method returns true only if segments have proper intersection. If set true, returns also true if segments have improper intersection (end of one segment on other segment ; continuous segments).
* \returns Whether the segments intersect
* * Example:
* \code{.py}
* ret = QgsGeometryUtils.segmentIntersection( QgsPoint( 0, 0 ), QgsPoint( 0, 1 ), QgsPoint( 1, 1 ), QgsPoint( 1, 0 ) )
* ret[0], ret[1].asWkt(), ret[2]
* # Whether the segments intersect, the intersection point, is intersect
* # (False, 'Point (0 0)', False)
* ret = QgsGeometryUtils.segmentIntersection( QgsPoint( 0, 0 ), QgsPoint( 0, 5 ), QgsPoint( 0, 5 ), QgsPoint( 1, 5 ) )
* ret[0], ret[1].asWkt(), ret[2]
* # (False, 'Point (0 5)', True)
* ret = QgsGeometryUtils.segmentIntersection( QgsPoint( 0, 0 ), QgsPoint( 0, 5 ), QgsPoint( 0, 5 ), QgsPoint( 1, 5 ), acceptImproperIntersection=True )
* ret[0], ret[1].asWkt(), ret[2]
* # (True, 'Point (0 5)', True)
* ret = QgsGeometryUtils.segmentIntersection( QgsPoint( 0, 0 ), QgsPoint( 0, 5 ), QgsPoint( 0, 2 ), QgsPoint( 1, 5 ) )
* ret[0], ret[1].asWkt(), ret[2]
* # (False, 'Point (0 2)', True)
* ret = QgsGeometryUtils.segmentIntersection( QgsPoint( 0, 0 ), QgsPoint( 0, 5 ), QgsPoint( 0, 2 ), QgsPoint( 1, 5 ), acceptImproperIntersection=True )
* ret[0], ret[1].asWkt(), ret[2]
* # (True, 'Point (0 2)', True)
* ret = QgsGeometryUtils.segmentIntersection( QgsPoint( 0, -5 ), QgsPoint( 0, 5 ), QgsPoint( 2, 0 ), QgsPoint( -1, 0 ) )
* ret[0], ret[1].asWkt(), ret[2]
* # (True, 'Point (0 0)', True)
* \endcode
*/
static bool segmentIntersection( const QgsPoint &p1, const QgsPoint &p2, const QgsPoint &q1, const QgsPoint &q2, QgsPoint &inter SIP_OUT, double tolerance );
static bool segmentIntersection( const QgsPoint &p1, const QgsPoint &p2, const QgsPoint &q1, const QgsPoint &q2, QgsPoint &intersectionPoint SIP_OUT, bool &isIntersection SIP_OUT, const double tolerance = 1e-8, bool acceptImproperIntersection = false );

/**
* \brief Project the point on a segment
@@ -495,14 +495,15 @@ QVector<QgsLineString> QgsTriangle::bisectors( double lengthTolerance ) const
QgsLineString bis3;
QgsPoint incenter = inscribedCenter();
QgsPoint out;
bool intersection = false;

QgsGeometryUtils::segmentIntersection( vertexAt( 0 ), incenter, vertexAt( 1 ), vertexAt( 2 ), out, lengthTolerance );
QgsGeometryUtils::segmentIntersection( vertexAt( 0 ), incenter, vertexAt( 1 ), vertexAt( 2 ), out, intersection, lengthTolerance );
bis1.setPoints( QgsPointSequence() << vertexAt( 0 ) << out );

QgsGeometryUtils::segmentIntersection( vertexAt( 1 ), incenter, vertexAt( 0 ), vertexAt( 2 ), out, lengthTolerance );
QgsGeometryUtils::segmentIntersection( vertexAt( 1 ), incenter, vertexAt( 0 ), vertexAt( 2 ), out, intersection, lengthTolerance );
bis2.setPoints( QgsPointSequence() << vertexAt( 1 ) << out );

QgsGeometryUtils::segmentIntersection( vertexAt( 2 ), incenter, vertexAt( 0 ), vertexAt( 1 ), out, lengthTolerance );
QgsGeometryUtils::segmentIntersection( vertexAt( 2 ), incenter, vertexAt( 0 ), vertexAt( 1 ), out, intersection, lengthTolerance );
bis3.setPoints( QgsPointSequence() << vertexAt( 2 ) << out );

bis.append( bis1 );
@@ -529,7 +530,8 @@ QgsPoint QgsTriangle::orthocenter( double lengthTolerance ) const
return QgsPoint();
QVector<QgsLineString> alt = altitudes();
QgsPoint ortho;
QgsGeometryUtils::segmentIntersection( alt.at( 0 ).pointN( 0 ), alt.at( 0 ).pointN( 1 ), alt.at( 1 ).pointN( 0 ), alt.at( 1 ).pointN( 1 ), ortho, lengthTolerance );
bool intersection;
QgsGeometryUtils::segmentIntersection( alt.at( 0 ).pointN( 0 ), alt.at( 0 ).pointN( 1 ), alt.at( 1 ).pointN( 0 ), alt.at( 1 ).pointN( 1 ), ortho, intersection, lengthTolerance );

return ortho;
}

0 comments on commit 4ff72de

Please sign in to comment.
You can’t perform that action at this time.