Skip to content
Permalink
Browse files

[FEATURE] New geometry API call to return a curve substring

Returns a new curve representing a substring of a curve, from
a start distance and end distance.

If z or m values are present, the output z and m will be interpolated using
the existing vertices' z or m values.

Handles curved geometries without loss.
  • Loading branch information
nyalldawson committed Aug 10, 2018
1 parent 5b0bdbd commit 513bcb68e48faff5546b7c7ed02e7f7cc68eb9b3
@@ -143,6 +143,8 @@ Sets the circular string's points

virtual QgsCircularString *reversed() const /Factory/;

virtual QgsCircularString *curveSubstring( double startDistance, double endDistance ) const /Factory/;

virtual bool addZValue( double zValue = 0 );

virtual bool addMValue( double mValue = 0 );
@@ -144,6 +144,8 @@ Appends first point if not already closed.

virtual QgsCompoundCurve *reversed() const /Factory/;

virtual QgsCompoundCurve *curveSubstring( double startDistance, double endDistance ) const /Factory/;


virtual bool addZValue( double zValue = 0 );

@@ -180,6 +180,20 @@ Returns the y-coordinate of the specified node in the line string.
virtual QPolygonF asQPolygonF() const;
%Docstring
Returns a QPolygonF representing the points.
%End

virtual QgsCurve *curveSubstring( double startDistance, double endDistance ) const = 0 /Factory/;
%Docstring
Returns a new curve representing a substring of this curve.

The ``startDistance`` and ``endDistance`` arguments specify the length along the curve
which the substring should start and end at. If the ``endDistance`` is greater than the
total length of the curve then any "extra" length will be ignored.

If z or m values are present, the output z and m will be interpolated using
the existing vertices' z or m values.

.. versionadded:: 3.4
%End

double straightDistance2d() const;
@@ -283,6 +283,8 @@ of the curve.

virtual QgsLineString *reversed() const /Factory/;

virtual QgsLineString *curveSubstring( double startDistance, double endDistance ) const /Factory/;


virtual double closestSegment( const QgsPoint &pt, QgsPoint &segmentPt /Out/, QgsVertexId &vertexAfter /Out/, int *leftOf /Out/ = 0, double epsilon = 4 * DBL_EPSILON ) const;

@@ -1178,6 +1178,133 @@ QgsCircularString *QgsCircularString::reversed() const
return copy;
}

QgsCircularString *QgsCircularString::curveSubstring( double startDistance, double endDistance ) const
{
if ( startDistance < 0 && endDistance < 0 )
return createEmptyWithSameType();

endDistance = std::max( startDistance, endDistance );

double distanceTraversed = 0;
const int totalPoints = numPoints();
if ( totalPoints == 0 )
return clone();

QVector< QgsPoint > substringPoints;

QgsWkbTypes::Type pointType = QgsWkbTypes::Point;
if ( is3D() )
pointType = QgsWkbTypes::PointZ;
if ( isMeasure() )
pointType = QgsWkbTypes::addM( pointType );

const double *x = mX.constData();
const double *y = mY.constData();
const double *z = is3D() ? mZ.constData() : nullptr;
const double *m = isMeasure() ? mM.constData() : nullptr;

double prevX = *x++;
double prevY = *y++;
double prevZ = z ? *z++ : 0.0;
double prevM = m ? *m++ : 0.0;
bool foundStart = false;

if ( qgsDoubleNear( startDistance, 0.0 ) || startDistance < 0 )
{
substringPoints << QgsPoint( pointType, prevX, prevY, prevZ, prevM );
foundStart = true;
}

substringPoints.reserve( totalPoints );

for ( int i = 0; i < ( totalPoints - 2 ) ; i += 2 )
{
double x1 = prevX;
double y1 = prevY;
double z1 = prevZ;
double m1 = prevM;

double x2 = *x++;
double y2 = *y++;
double z2 = z ? *z++ : 0.0;
double m2 = m ? *m++ : 0.0;

double x3 = *x++;
double y3 = *y++;
double z3 = z ? *z++ : 0.0;
double m3 = m ? *m++ : 0.0;

bool addedSegmentEnd = false;
const double segmentLength = QgsGeometryUtils::circleLength( x1, y1, x2, y2, x3, y3 );
if ( distanceTraversed < startDistance && distanceTraversed + segmentLength > startDistance )
{
// start point falls on this segment
const double distanceToStart = startDistance - distanceTraversed;
const QgsPoint startPoint = QgsGeometryUtils::interpolatePointOnArc( QgsPoint( pointType, x1, y1, z1, m1 ),
QgsPoint( pointType, x2, y2, z2, m2 ),
QgsPoint( pointType, x3, y3, z3, m3 ), distanceToStart );

// does end point also fall on this segment?
const bool endPointOnSegment = distanceTraversed + segmentLength > endDistance;
if ( endPointOnSegment )
{
const double distanceToEnd = endDistance - distanceTraversed;
const double midPointDistance = ( distanceToEnd - distanceToStart ) * 0.5 + distanceToStart;
substringPoints << startPoint
<< QgsGeometryUtils::interpolatePointOnArc( QgsPoint( pointType, x1, y1, z1, m1 ),
QgsPoint( pointType, x2, y2, z2, m2 ),
QgsPoint( pointType, x3, y3, z3, m3 ), midPointDistance )
<< QgsGeometryUtils::interpolatePointOnArc( QgsPoint( pointType, x1, y1, z1, m1 ),
QgsPoint( pointType, x2, y2, z2, m2 ),
QgsPoint( pointType, x3, y3, z3, m3 ), distanceToEnd );
addedSegmentEnd = true;
}
else
{
const double midPointDistance = ( segmentLength - distanceToStart ) * 0.5 + distanceToStart;
substringPoints << startPoint
<< QgsGeometryUtils::interpolatePointOnArc( QgsPoint( pointType, x1, y1, z1, m1 ),
QgsPoint( pointType, x2, y2, z2, m2 ),
QgsPoint( pointType, x3, y3, z3, m3 ), midPointDistance )
<< QgsPoint( pointType, x3, y3, z3, m3 );
addedSegmentEnd = true;
}
foundStart = true;
}
if ( !addedSegmentEnd && foundStart && ( distanceTraversed + segmentLength > endDistance ) )
{
// end point falls on this segment
const double distanceToEnd = endDistance - distanceTraversed;
// add mid point, at half way along this arc, then add the interpolated end point
substringPoints << QgsGeometryUtils::interpolatePointOnArc( QgsPoint( pointType, x1, y1, z1, m1 ),
QgsPoint( pointType, x2, y2, z2, m2 ),
QgsPoint( pointType, x3, y3, z3, m3 ), distanceToEnd / 2.0 )

<< QgsGeometryUtils::interpolatePointOnArc( QgsPoint( pointType, x1, y1, z1, m1 ),
QgsPoint( pointType, x2, y2, z2, m2 ),
QgsPoint( pointType, x3, y3, z3, m3 ), distanceToEnd );
}
else if ( !addedSegmentEnd && foundStart )
{
substringPoints << QgsPoint( pointType, x2, y2, z2, m2 )
<< QgsPoint( pointType, x3, y3, z3, m3 );
}

distanceTraversed += segmentLength;
if ( distanceTraversed > endDistance )
break;

prevX = x3;
prevY = y3;
prevZ = z3;
prevM = m3;
}

std::unique_ptr< QgsCircularString > result = qgis::make_unique< QgsCircularString >();
result->setPoints( substringPoints );
return result.release();
}

bool QgsCircularString::addZValue( double zValue )
{
if ( QgsWkbTypes::hasZ( mWkbType ) )
@@ -118,6 +118,7 @@ class CORE_EXPORT QgsCircularString: public QgsCurve
double vertexAngle( QgsVertexId vertex ) const override;
double segmentLength( QgsVertexId startVertex ) const override;
QgsCircularString *reversed() const override SIP_FACTORY;
QgsCircularString *curveSubstring( double startDistance, double endDistance ) const override SIP_FACTORY;
bool addZValue( double zValue = 0 ) override;
bool addMValue( double mValue = 0 ) override;
bool dropZValue() override;
@@ -864,6 +864,37 @@ QgsCompoundCurve *QgsCompoundCurve::reversed() const
return clone;
}

QgsCompoundCurve *QgsCompoundCurve::curveSubstring( double startDistance, double endDistance ) const
{
if ( startDistance < 0 && endDistance < 0 )
return createEmptyWithSameType();

endDistance = std::max( startDistance, endDistance );
std::unique_ptr< QgsCompoundCurve > substring = qgis::make_unique< QgsCompoundCurve >();

double distanceTraversed = 0;
for ( const QgsCurve *curve : mCurves )
{
const double thisCurveLength = curve->length();
if ( distanceTraversed + thisCurveLength < startDistance )
{
// keep going - haven't found start yet, so no need to include this curve at all
}
else
{
std::unique_ptr< QgsCurve > part( curve->curveSubstring( startDistance - distanceTraversed, endDistance - distanceTraversed ) );
if ( part )
substring->addCurve( part.release() );
}

distanceTraversed += thisCurveLength;
if ( distanceTraversed > endDistance )
break;
}

return substring.release();
}

bool QgsCompoundCurve::addZValue( double zValue )
{
if ( QgsWkbTypes::hasZ( mWkbType ) )
@@ -115,6 +115,7 @@ class CORE_EXPORT QgsCompoundCurve: public QgsCurve
double vertexAngle( QgsVertexId vertex ) const override;
double segmentLength( QgsVertexId startVertex ) const override;
QgsCompoundCurve *reversed() const override SIP_FACTORY;
QgsCompoundCurve *curveSubstring( double startDistance, double endDistance ) const override SIP_FACTORY;

bool addZValue( double zValue = 0 ) override;
bool addMValue( double mValue = 0 ) override;
@@ -166,6 +166,20 @@ class CORE_EXPORT QgsCurve: public QgsAbstractGeometry
*/
virtual QPolygonF asQPolygonF() const;

/**
* Returns a new curve representing a substring of this curve.
*
* The \a startDistance and \a endDistance arguments specify the length along the curve
* which the substring should start and end at. If the \a endDistance is greater than the
* total length of the curve then any "extra" length will be ignored.
*
* If z or m values are present, the output z and m will be interpolated using
* the existing vertices' z or m values.
*
* \since QGIS 3.4
*/
virtual QgsCurve *curveSubstring( double startDistance, double endDistance ) const = 0 SIP_FACTORY;

/**
* Returns the straight distance of the curve, i.e. the direct/euclidean distance
* between the first and last vertex of the curve. (Also known as
@@ -711,6 +711,96 @@ QgsLineString *QgsLineString::reversed() const
return copy;
}

QgsLineString *QgsLineString::curveSubstring( double startDistance, double endDistance ) const
{
if ( startDistance < 0 && endDistance < 0 )
return createEmptyWithSameType();

endDistance = std::max( startDistance, endDistance );

double distanceTraversed = 0;
const int totalPoints = numPoints();
if ( totalPoints == 0 )
return clone();

QVector< QgsPoint > substringPoints;

QgsWkbTypes::Type pointType = QgsWkbTypes::Point;
if ( is3D() )
pointType = QgsWkbTypes::PointZ;
if ( isMeasure() )
pointType = QgsWkbTypes::addM( pointType );

const double *x = mX.constData();
const double *y = mY.constData();
const double *z = is3D() ? mZ.constData() : nullptr;
const double *m = isMeasure() ? mM.constData() : nullptr;

double prevX = *x++;
double prevY = *y++;
double prevZ = z ? *z++ : 0.0;
double prevM = m ? *m++ : 0.0;
bool foundStart = false;

if ( qgsDoubleNear( startDistance, 0.0 ) || startDistance < 0 )
{
substringPoints << QgsPoint( pointType, prevX, prevY, prevZ, prevM );
foundStart = true;
}

substringPoints.reserve( totalPoints );

for ( int i = 1; i < totalPoints; ++i )
{
double thisX = *x++;
double thisY = *y++;
double thisZ = z ? *z++ : 0.0;
double thisM = m ? *m++ : 0.0;

const double segmentLength = std::sqrt( ( thisX - prevX ) * ( thisX - prevX ) + ( thisY - prevY ) * ( thisY - prevY ) );
if ( distanceTraversed < startDistance && distanceTraversed + segmentLength > startDistance )
{
// start point falls on this segment
const double distanceToStart = startDistance - distanceTraversed;
double startX, startY;
double startZ = 0;
double startM = 0;
QgsGeometryUtils::pointOnLineWithDistance( prevX, prevY, thisX, thisY, distanceToStart, startX, startY,
z ? &prevZ : nullptr, z ? &thisZ : nullptr, z ? &startZ : nullptr,
m ? &prevM : nullptr, m ? &thisM : nullptr, m ? &startM : nullptr );
substringPoints << QgsPoint( pointType, startX, startY, startZ, startM );
foundStart = true;
}
if ( foundStart && ( distanceTraversed + segmentLength > endDistance ) )
{
// end point falls on this segment
const double distanceToEnd = endDistance - distanceTraversed;
double endX, endY;
double endZ = 0;
double endM = 0;
QgsGeometryUtils::pointOnLineWithDistance( prevX, prevY, thisX, thisY, distanceToEnd, endX, endY,
z ? &prevZ : nullptr, z ? &thisZ : nullptr, z ? &endZ : nullptr,
m ? &prevM : nullptr, m ? &thisM : nullptr, m ? &endM : nullptr );
substringPoints << QgsPoint( pointType, endX, endY, endZ, endM );
}
else if ( foundStart )
{
substringPoints << QgsPoint( pointType, thisX, thisY, thisZ, thisM );
}

distanceTraversed += segmentLength;
if ( distanceTraversed > endDistance )
break;

prevX = thisX;
prevY = thisY;
prevZ = thisZ;
prevM = thisM;
}

return new QgsLineString( substringPoints );
}

/***************************************************************************
* This class is considered CRITICAL and any change MUST be accompanied with
* full unit tests.
@@ -308,6 +308,7 @@ class CORE_EXPORT QgsLineString: public QgsCurve
bool deleteVertex( QgsVertexId position ) override;

QgsLineString *reversed() const override SIP_FACTORY;
QgsLineString *curveSubstring( double startDistance, double endDistance ) const override SIP_FACTORY;

double closestSegment( const QgsPoint &pt, QgsPoint &segmentPt SIP_OUT, QgsVertexId &vertexAfter SIP_OUT, int *leftOf SIP_OUT = nullptr, double epsilon = 4 * std::numeric_limits<double>::epsilon() ) const override;
bool pointAt( int node, QgsPoint &point, QgsVertexId::VertexType &type ) const override;

0 comments on commit 513bcb6

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