From f439a649d64db181b32f2dbca05f8b1af6158af8 Mon Sep 17 00:00:00 2001 From: Nyall Dawson Date: Wed, 25 May 2022 15:20:35 +1000 Subject: [PATCH] Optimise 3d clipping to work directly with coordinate arrays, so that we avoid the (signficant) overhead of point conversion and sequence allocations --- .../geometry/qgslinestring.sip.in | 5 + python/core/auto_generated/qgsclipper.sip.in | 20 -- src/core/geometry/qgslinestring.h | 41 ++++ src/core/qgsclipper.cpp | 115 ++++++--- src/core/qgsclipper.h | 169 +++++++------ src/core/symbology/qgssymbol.cpp | 223 ++++++++++++++---- tests/src/core/testqgsclipper.cpp | 60 ++--- 7 files changed, 430 insertions(+), 203 deletions(-) diff --git a/python/core/auto_generated/geometry/qgslinestring.sip.in b/python/core/auto_generated/geometry/qgslinestring.sip.in index be905bff830..152af6671c4 100644 --- a/python/core/auto_generated/geometry/qgslinestring.sip.in +++ b/python/core/auto_generated/geometry/qgslinestring.sip.in @@ -358,6 +358,11 @@ corresponds to the last point in the line. + + + + + double zAt( int index ) const; %Docstring Returns the z-coordinate of the specified node in the line string. diff --git a/python/core/auto_generated/qgsclipper.sip.in b/python/core/auto_generated/qgsclipper.sip.in index 044cb076252..88d77e16b8a 100644 --- a/python/core/auto_generated/qgsclipper.sip.in +++ b/python/core/auto_generated/qgsclipper.sip.in @@ -77,27 +77,7 @@ Trims the given polygon to a rectangular box, by modifying the given polygon in :param clipRect: clipping rectangle %End - static void trimPolygon( QgsLineString &polygon, const QgsBox3d &clipRect ); -%Docstring -Trims the given polygon to a pseudo 3D box, by modifying the given polygon in place. -:param polygon: polygon as 2D or 3D coordinates -:param clipRect: clipping 3D box - -.. versionadded:: 3.22 -%End - - static QgsLineString clipped3dLine( const QgsCurve &curve, const QgsBox3d &clipExtent ); -%Docstring -Takes a ``curve`` with 3D coordinates and clips it to clipExtent - -:param curve: the linestring to clip -:param clipExtent: clipping bounds - -:return: clipped line coordinates - -.. versionadded:: 3.22 -%End static QPolygonF clippedLine( const QgsCurve &curve, const QgsRectangle &clipExtent ); %Docstring diff --git a/src/core/geometry/qgslinestring.h b/src/core/geometry/qgslinestring.h index 3e7977dccf8..67c2808e57a 100644 --- a/src/core/geometry/qgslinestring.h +++ b/src/core/geometry/qgslinestring.h @@ -458,6 +458,47 @@ class CORE_EXPORT QgsLineString: public QgsCurve return mM.constData(); } + /** + * Returns the x vertex values as a vector. + * \note Not available in Python bindings + * \since QGIS 3.26 + */ + QVector< double > xVector() const SIP_SKIP + { + return mX; + } + + /** + * Returns the y vertex values as a vector. + * \note Not available in Python bindings + * \since QGIS 3.26 + */ + QVector< double > yVector() const SIP_SKIP + { + return mY; + } + + /** + * Returns the z vertex values as a vector. + * \note Not available in Python bindings + * \since QGIS 3.26 + */ + QVector< double > zVector() const SIP_SKIP + { + return mZ; + } + + /** + * Returns the m vertex values as a vector. + * \note Not available in Python bindings + * \since QGIS 3.26 + */ + QVector< double > mVector() const SIP_SKIP + { + return mM; + } + + #ifndef SIP_RUN /** diff --git a/src/core/qgsclipper.cpp b/src/core/qgsclipper.cpp index 6830e74ebe4..c0af111393b 100644 --- a/src/core/qgsclipper.cpp +++ b/src/core/qgsclipper.cpp @@ -37,22 +37,29 @@ const double QgsClipper::MIN_Y = -16000; const double QgsClipper::SMALL_NUM = 1e-12; -QgsLineString QgsClipper::clipped3dLine( const QgsCurve &curve, const QgsBox3d &clipExtent ) +void QgsClipper::clipped3dLine( const QVector< double > &xIn, const QVector< double > &yIn, const QVector &zIn, QVector &x, QVector &y, QVector &z, const QgsBox3d &clipExtent ) { double p0x, p0y, p0z, p1x = 0.0, p1y = 0.0, p1z = 0.0; //original coordinates double p1x_c, p1y_c, p1z_c; //clipped end coordinates double lastClipX = 0.0, lastClipY = 0.0, lastClipZ = 0.0; // last successfully clipped coordinates - QgsPointSequence seq; + const int nPoints = xIn.size(); - for ( QgsAbstractGeometry::vertex_iterator ite = curve.vertices_begin() ; ite != curve.vertices_end(); ++ite ) + x.reserve( nPoints ); + y.reserve( nPoints ); + z.reserve( nPoints ); + + const double *sourceX = xIn.data(); + const double *sourceY = yIn.data(); + const double *sourceZ = zIn.data(); + + for ( int i = 0; i < nPoints; ++i ) { - QgsPoint curPoint = *ite; - if ( ite == curve.vertices_begin() ) + if ( i == 0 ) { - p1x = curPoint.x(); - p1y = curPoint.y(); - p1z = curPoint.z(); + p1x = *sourceX++; + p1y = *sourceY++; + p1z = *sourceZ++; } else { @@ -60,9 +67,9 @@ QgsLineString QgsClipper::clipped3dLine( const QgsCurve &curve, const QgsBox3d & p0y = p1y; p0z = p1z; - p1x = curPoint.x(); - p1y = curPoint.y(); - p1z = curPoint.z(); + p1x = *sourceX++; + p1y = *sourceY++; + p1z = *sourceZ++; p1x_c = p1x; p1y_c = p1y; @@ -71,31 +78,33 @@ QgsLineString QgsClipper::clipped3dLine( const QgsCurve &curve, const QgsBox3d & // TODO: should be in 3D if ( clipLineSegment( clipExtent, p0x, p0y, p0z, p1x_c, p1y_c, p1z_c ) ) { - bool newLine = !seq.isEmpty() && ( !qgsDoubleNear( p0x, lastClipX ) - || !qgsDoubleNear( p0y, lastClipY ) - || !qgsDoubleNear( p0z, lastClipZ ) ); + bool newLine = !x.isEmpty() && ( !qgsDoubleNear( p0x, lastClipX ) + || !qgsDoubleNear( p0y, lastClipY ) + || !qgsDoubleNear( p0z, lastClipZ ) ); if ( newLine ) { //add edge points to connect old and new line // TODO: should be (really) in 3D - connectSeparatedLines( lastClipX, lastClipY, lastClipZ, p0x, p0y, p0z, clipExtent, seq ); + connectSeparatedLines( lastClipX, lastClipY, lastClipZ, p0x, p0y, p0z, clipExtent, x, y, z ); } - if ( seq.isEmpty() || newLine ) + if ( x.isEmpty() || newLine ) { //add first point - seq << QgsPoint( p0x, p0y, p0z ) ; + x << p0x; + y << p0y; + z << p0z; } //add second point lastClipX = p1x_c; lastClipY = p1y_c; lastClipZ = p1z_c; - seq << QgsPoint( p1x_c, p1y_c, p1z_c ) ; + x << p1x_c; + y << p1y_c; + z << p1z_c; } } } - - return QgsLineString( seq ); } QPolygonF QgsClipper::clippedLine( const QgsCurve &curve, const QgsRectangle &clipExtent ) @@ -258,7 +267,7 @@ void QgsClipper::connectSeparatedLines( double x0, double y0, double x1, double } void QgsClipper::connectSeparatedLines( double x0, double y0, double z0, double x1, double y1, double z1, - const QgsBox3d &clipRect, QgsPointSequence &pts ) + const QgsBox3d &clipRect, QVector< double > &ptsX, QVector< double > &ptsY, QVector &ptsZ ) { // TODO: really relevant and sufficient? double meanZ = ( z0 + z1 ) / 2.0; @@ -272,18 +281,26 @@ void QgsClipper::connectSeparatedLines( double x0, double y0, double z0, double } else if ( qgsDoubleNear( y1, clipRect.yMaximum() ) ) { - pts << QgsPoint( clipRect.xMinimum(), clipRect.yMaximum(), meanZ ); + ptsX << clipRect.xMinimum(); + ptsY << clipRect.yMaximum(); + ptsZ << meanZ; return; } else if ( qgsDoubleNear( x1, clipRect.xMaximum() ) ) { - pts << QgsPoint( clipRect.xMinimum(), clipRect.yMaximum(), meanZ ); - pts << QgsPoint( clipRect.xMaximum(), clipRect.yMaximum(), meanZ ); + ptsX << clipRect.xMinimum(); + ptsY << clipRect.yMaximum(); + ptsZ << meanZ; + ptsX << clipRect.xMaximum(); + ptsY << clipRect.yMaximum(); + ptsZ << meanZ; return; } else if ( qgsDoubleNear( y1, clipRect.yMinimum() ) ) { - pts << QgsPoint( clipRect.xMinimum(), clipRect.yMinimum(), meanZ ); + ptsX << clipRect.xMinimum(); + ptsY << clipRect.yMinimum(); + ptsZ << meanZ; return; } } @@ -295,18 +312,26 @@ void QgsClipper::connectSeparatedLines( double x0, double y0, double z0, double } else if ( qgsDoubleNear( x1, clipRect.xMaximum() ) ) { - pts << QgsPoint( clipRect.xMaximum(), clipRect.yMaximum(), meanZ ); + ptsX << clipRect.xMaximum(); + ptsY << clipRect.yMaximum(); + ptsZ << meanZ; return; } else if ( qgsDoubleNear( y1, clipRect.yMinimum() ) ) { - pts << QgsPoint( clipRect.xMaximum(), clipRect.yMaximum(), meanZ ); - pts << QgsPoint( clipRect.xMaximum(), clipRect.yMinimum(), meanZ ); + ptsX << clipRect.xMaximum(); + ptsY << clipRect.yMaximum(); + ptsZ << meanZ; + ptsX << clipRect.xMaximum(); + ptsY << clipRect.yMinimum(); + ptsZ << meanZ; return; } else if ( qgsDoubleNear( x1, clipRect.xMinimum() ) ) { - pts << QgsPoint( clipRect.xMinimum(), clipRect.yMaximum(), meanZ ); + ptsX << clipRect.xMinimum(); + ptsY << clipRect.yMaximum(); + ptsZ << meanZ; return; } } @@ -318,18 +343,26 @@ void QgsClipper::connectSeparatedLines( double x0, double y0, double z0, double } else if ( qgsDoubleNear( y1, clipRect.yMinimum() ) ) { - pts << QgsPoint( clipRect.xMaximum(), clipRect.yMinimum(), meanZ ); + ptsX << clipRect.xMaximum(); + ptsY << clipRect.yMinimum(); + ptsZ << meanZ; return; } else if ( qgsDoubleNear( x1, clipRect.xMinimum() ) ) { - pts << QgsPoint( clipRect.xMaximum(), clipRect.yMinimum(), meanZ ); - pts << QgsPoint( clipRect.xMinimum(), clipRect.yMinimum(), meanZ ); + ptsX << clipRect.xMaximum(); + ptsY << clipRect.yMinimum(); + ptsZ << meanZ; + ptsX << clipRect.xMinimum(); + ptsY << clipRect.yMinimum(); + ptsZ << meanZ; return; } else if ( qgsDoubleNear( y1, clipRect.yMaximum() ) ) { - pts << QgsPoint( clipRect.xMaximum(), clipRect.yMaximum(), meanZ ); + ptsX << clipRect.xMaximum(); + ptsY << clipRect.yMaximum(); + ptsZ << meanZ; return; } } @@ -341,18 +374,26 @@ void QgsClipper::connectSeparatedLines( double x0, double y0, double z0, double } else if ( qgsDoubleNear( x1, clipRect.xMinimum() ) ) { - pts << QgsPoint( clipRect.xMinimum(), clipRect.yMinimum(), meanZ ); + ptsX << clipRect.xMinimum(); + ptsY << clipRect.yMinimum(); + ptsZ << meanZ; return; } else if ( qgsDoubleNear( y1, clipRect.yMaximum() ) ) { - pts << QgsPoint( clipRect.xMinimum(), clipRect.yMinimum(), meanZ ); - pts << QgsPoint( clipRect.xMinimum(), clipRect.yMaximum(), meanZ ); + ptsX << clipRect.xMinimum(); + ptsY << clipRect.yMinimum(); + ptsZ << meanZ; + ptsX << clipRect.xMinimum(); + ptsY << clipRect.yMaximum(); + ptsZ << meanZ; return; } else if ( qgsDoubleNear( x1, clipRect.xMaximum() ) ) { - pts << QgsPoint( clipRect.xMaximum(), clipRect.yMinimum(), meanZ ); + ptsX << clipRect.xMaximum(); + ptsY << clipRect.yMinimum(); + ptsZ << meanZ; return; } } diff --git a/src/core/qgsclipper.h b/src/core/qgsclipper.h index 2e58d623022..f3af128e72e 100644 --- a/src/core/qgsclipper.h +++ b/src/core/qgsclipper.h @@ -85,8 +85,8 @@ class CORE_EXPORT QgsClipper XMin, YMax, YMin, - ZMax, //!< Maximum Z (since QGIS 3.22) - ZMin, //!< Minimum Z (since QGIS 3.22) + ZMax, //!< Maximum Z (since QGIS 3.26) + ZMin, //!< Minimum Z (since QGIS 3.26) }; SIP_IF_FEATURE( !ARM ) // Not available on ARM sip bindings because of qreal issues @@ -114,22 +114,19 @@ class CORE_EXPORT QgsClipper static void trimPolygon( QPolygonF &pts, const QgsRectangle &clipRect ); /** - * Trims the given polygon to a pseudo 3D box, by modifying the given polygon in place. + * Trims a polygon consisting of the specified \a x, \a y and \a z values to a pseudo 3D box, by modifying the arrays in place. * - * \param polygon polygon as 2D or 3D coordinates - * \param clipRect clipping 3D box - * \since QGIS 3.22 + * \note Not available in Python bindings + * \since QGIS 3.26 */ - static void trimPolygon( QgsLineString &polygon, const QgsBox3d &clipRect ); + static void trimPolygon( QVector< double > &x, QVector< double > &y, QVector< double> &z, const QgsBox3d &clipRect ) SIP_SKIP; /** - * Takes a \a curve with 3D coordinates and clips it to clipExtent - * \param curve the linestring to clip - * \param clipExtent clipping bounds - * \returns clipped line coordinates - * \since QGIS 3.22 + * Takes a line with 3D coordinates and clips it to clipExtent. + * \note Not available in Python bindings + * \since QGIS 3.26 */ - static QgsLineString clipped3dLine( const QgsCurve &curve, const QgsBox3d &clipExtent ); + static void clipped3dLine( const QVector< double > &xIn, const QVector< double > &yIn, const QVector &zIn, QVector< double > &x, QVector< double > &y, QVector< double > &z, const QgsBox3d &clipExtent ) SIP_SKIP; /** * Takes a linestring and clips it to clipExtent @@ -183,14 +180,15 @@ class CORE_EXPORT QgsClipper static inline void trimPolygonToBoundary( const QPolygonF &inPts, QPolygonF &outPts, const QgsRectangle &rect, Boundary b, double boundaryValue ); - static inline void trimPolygonToBoundary( const QgsPointSequence &inPts, QgsPointSequence &outPts, const QgsBox3d &rect, Boundary b, double boundaryValue ); + static inline void trimPolygonToBoundary( const QVector &inX, const QVector &inY, const QVector< double > &inZ, + QVector &outX, QVector &outY, QVector &outZ, const QgsBox3d &rect, Boundary boundary, double boundaryValue ); // Determines if a point is inside or outside the given boundary static inline bool inside( double x, double y, Boundary b ); static inline bool inside( QPointF pt, Boundary b, double val ); - static inline bool inside( QgsPoint pt, Boundary b, double val ); + static inline bool inside( double x, double y, double z, Boundary boundary, double boundaryValue ); // Calculates the intersection point between a line defined by a // (x1, y1), and (x2, y2) and the given boundary @@ -202,9 +200,9 @@ class CORE_EXPORT QgsClipper QPointF pt2, Boundary b, const QgsRectangle &rect ); - static inline QgsPoint intersectRect( QgsPoint pt1, - QgsPoint pt2, - Boundary b, const QgsBox3d &rect ); + static inline void intersectRect( double x1, double y1, double z1, double x2, double y2, double z2, + double &xOut, double &yOut, double &zOut, + Boundary boundary, const QgsBox3d &rect ); static bool clipLineSegment( const QgsBox3d &extent, double &x0, double &y0, double &z0, double &x1, double &y1, double &z1 ); @@ -229,10 +227,12 @@ class CORE_EXPORT QgsClipper * \param y1 y-coordinate of the second line start * \param z1 z-coordinate of the second line start * \param clipRect clip box - * \param pts: in/out array of clipped points + * \param ptsX: in/out array of clipped points x + * \param ptsY: in/out array of clipped points y + * \param ptsZ: in/out array of clipped points z */ static void connectSeparatedLines( double x0, double y0, double z0, double x1, double y1, double z1, - const QgsBox3d &clipRect, QgsPointSequence &pts ); + const QgsBox3d &clipRect, QVector< double > &ptsX, QVector< double > &ptsY, QVector &ptsZ ); //low level clip methods for fast clip algorithm static void clipStartTop( double &x0, double &y0, double x1, double y1, double yMax ); @@ -296,27 +296,41 @@ inline void QgsClipper::trimPolygon( QPolygonF &pts, const QgsRectangle &clipRec trimPolygonToBoundary( tmpPts, pts, clipRect, YMin, clipRect.yMinimum() ); } -inline void QgsClipper::trimPolygon( QgsLineString &pts, const QgsBox3d &clipRect ) +inline void QgsClipper::trimPolygon( QVector &x, QVector &y, QVector &z, const QgsBox3d &clipRect ) { - QgsPointSequence seqPts, seqTmpPts; - pts.points( seqPts ); + QVector< double > tempX; + QVector< double > tempY; + QVector< double > tempZ; + const int size = x.size(); + tempX.reserve( size ); + tempY.reserve( size ); + tempZ.reserve( size ); + + trimPolygonToBoundary( x, y, z, tempX, tempY, tempZ, clipRect, XMax, clipRect.xMaximum() ); + x.resize( 0 ); + y.resize( 0 ); + z.resize( 0 ); + trimPolygonToBoundary( tempX, tempY, tempZ, x, y, z, clipRect, YMax, clipRect.yMaximum() ); + tempX.resize( 0 ); + tempY.resize( 0 ); + tempZ.resize( 0 ); + trimPolygonToBoundary( x, y, z, tempX, tempY, tempZ, clipRect, XMin, clipRect.xMinimum() ); + x.resize( 0 ); + y.resize( 0 ); + z.resize( 0 ); + trimPolygonToBoundary( tempX, tempY, tempZ, x, y, z, clipRect, YMin, clipRect.yMinimum() ); - trimPolygonToBoundary( seqPts, seqTmpPts, clipRect, XMax, clipRect.xMaximum() ); - seqPts.clear(); - trimPolygonToBoundary( seqTmpPts, seqPts, clipRect, YMax, clipRect.yMaximum() ); - seqTmpPts.clear(); - trimPolygonToBoundary( seqPts, seqTmpPts, clipRect, XMin, clipRect.xMinimum() ); - seqPts.clear(); - trimPolygonToBoundary( seqTmpPts, seqPts, clipRect, YMin, clipRect.yMinimum() ); if ( !clipRect.is2d() ) { - seqTmpPts.clear(); - trimPolygonToBoundary( seqPts, seqTmpPts, clipRect, ZMax, clipRect.zMaximum() ); - seqPts.clear(); - trimPolygonToBoundary( seqTmpPts, seqPts, clipRect, ZMin, clipRect.zMinimum() ); + tempX.resize( 0 ); + tempY.resize( 0 ); + tempZ.resize( 0 ); + trimPolygonToBoundary( x, y, z, tempX, tempY, tempZ, clipRect, ZMax, clipRect.zMaximum() ); + x.resize( 0 ); + y.resize( 0 ); + z.resize( 0 ); + trimPolygonToBoundary( tempX, tempY, tempZ, x, y, z, clipRect, ZMin, clipRect.zMinimum() ); } - - pts.setPoints( seqPts ); } // An auxiliary function that is part of the polygon trimming @@ -424,40 +438,61 @@ inline void QgsClipper::trimPolygonToBoundary( const QPolygonF &inPts, QPolygonF } } -inline void QgsClipper::trimPolygonToBoundary( const QgsPointSequence &inPts, QgsPointSequence &outPts, const QgsBox3d &rect, Boundary b, double boundaryValue ) +inline void QgsClipper::trimPolygonToBoundary( const QVector &inX, const QVector &inY, const QVector< double > &inZ, + QVector &outX, QVector &outY, QVector &outZ, const QgsBox3d &rect, Boundary boundary, double boundaryValue ) { - QgsPoint inI1, inI2; + const double len = inX.length(); + double inI1X = inX.at( len - 1 ); // start with last point + double inI1Y = inY.at( len - 1 ); + double inI1Z = inZ.at( len - 1 ); - inI1 = inPts.at( inPts.length() - 1 ); // start with last point + double xTemp; + double yTemp; + double zTemp; // and compare to the first point initially. - for ( int i2 = 0; i2 < inPts.length() ; ++i2 ) + for ( int i2 = 0; i2 < len; ++i2 ) { - inI2 = inPts.at( i2 ); + double inI2X = inX.at( i2 ); + double inI2Y = inY.at( i2 ); + double inI2Z = inZ.at( i2 ); + // look at each edge of the polygon in turn - if ( inside( inI2, b, boundaryValue ) ) // end point of edge is inside boundary + if ( inside( inI2X, inI2Y, inI2Z, boundary, boundaryValue ) ) // end point of edge is inside boundary { - if ( inside( inI1, b, boundaryValue ) ) + if ( inside( inI1X, inI1Y, inI1Z, boundary, boundaryValue ) ) { - outPts << inI2 ; + outX << inI2X; + outY << inI2Y; + outZ << inI2Z; } else { // edge crosses into the boundary, so trim back to the boundary, and // store both ends of the new edge - outPts << intersectRect( inI1, inI2, b, rect ) ; - outPts << inI2 ; + intersectRect( inI1X, inI1Y, inI1Z, inI2X, inI2Y, inI2Z, xTemp, yTemp, zTemp, boundary, rect ) ; + outX << xTemp; + outY << yTemp; + outZ << zTemp; + outX << inI2X; + outY << inI2Y; + outZ << inI2Z; } } else // end point of edge is outside boundary { // start point is in boundary, so need to trim back - if ( inside( inI1, b, boundaryValue ) ) + if ( inside( inI1X, inI1Y, inI1Z, boundary, boundaryValue ) ) { - outPts << intersectRect( inI1, inI2, b, rect ) ; + intersectRect( inI1X, inI1Y, inI1Z, inI2X, inI2Y, inI2Z, xTemp, yTemp, zTemp, boundary, rect ); + outX << xTemp; + outY << yTemp; + outZ << zTemp; } } - inI1 = inPts.at( i2 ); + inI1X = inI2X; + inI1Y = inI2Y; + inI1Z = inI2Z; } } @@ -510,22 +545,22 @@ inline bool QgsClipper::inside( QPointF pt, Boundary b, double val ) return false; } -inline bool QgsClipper::inside( QgsPoint pt, Boundary b, double val ) +inline bool QgsClipper::inside( double x, double y, double z, Boundary boundary, double boundaryValue ) { - switch ( b ) + switch ( boundary ) { case XMax: // x < MAX_X is inside - return ( pt.x() < val ); + return ( x < boundaryValue ); case XMin: // x > MIN_X is inside - return ( pt.x() > val ); + return ( x > boundaryValue ); case YMax: // y < MAX_Y is inside - return ( pt.y() < val ); + return ( y < boundaryValue ); case YMin: // y > MIN_Y is inside - return ( pt.y() > val ); + return ( y > boundaryValue ); case ZMax: // z < MAX_Z is inside - return ( pt.z() < val ); + return ( z < boundaryValue ); case ZMin: // z > MIN_Z is inside - return ( pt.z() > val ); + return ( z > boundaryValue ); } return false; } @@ -629,20 +664,18 @@ inline QPointF QgsClipper::intersectRect( QPointF pt1, return QPointF( x1 + r * ( x2 - x1 ), y1 + r * ( y2 - y1 ) ); } -inline QgsPoint QgsClipper::intersectRect( QgsPoint pt1, - QgsPoint pt2, - Boundary b, const QgsBox3d &rect ) +inline void QgsClipper::intersectRect( const double x1, const double y1, const double z1, + const double x2, const double y2, const double z2, + double &xOut, double &yOut, double &zOut, + Boundary boundary, const QgsBox3d &rect ) { - // This function assumes that the two given points (x1, y1), and - // (x2, y2) cross the given boundary. Making this assumption allows + // This function assumes that the two given points (x1, y1, z1), and + // (x2, y2, z2) cross the given boundary. Making this assumption allows // some optimisations. double r_n = SMALL_NUM, r_d = SMALL_NUM; - const double x1 = pt1.x(), x2 = pt2.x(); - const double y1 = pt1.y(), y2 = pt2.y(); - const double z1 = pt1.z(), z2 = pt2.z(); - switch ( b ) + switch ( boundary ) { case XMax: // x = MAX_X boundary r_n = -( x1 - rect.xMaximum() ) * ( rect.yMaximum() - rect.yMinimum() ); @@ -675,7 +708,9 @@ inline QgsPoint QgsClipper::intersectRect( QgsPoint pt1, { r = r_n / r_d; } - return QgsPoint( x1 + r * ( x2 - x1 ), y1 + r * ( y2 - y1 ), z1 + r * ( z2 - z1 ) ); + xOut = x1 + r * ( x2 - x1 ); + yOut = y1 + r * ( y2 - y1 ); + zOut = z1 + r * ( z2 - z1 ); } inline void QgsClipper::clipStartTop( double &x0, double &y0, double x1, double y1, double yMax ) diff --git a/src/core/symbology/qgssymbol.cpp b/src/core/symbology/qgssymbol.cpp index c481b1fc1ed..85b671febf5 100644 --- a/src/core/symbology/qgssymbol.cpp +++ b/src/core/symbology/qgssymbol.cpp @@ -95,7 +95,9 @@ QPolygonF QgsSymbol::_getLineString3d( QgsRenderContext &context, const QgsCurve QgsCoordinateTransform ct = context.coordinateTransform(); const QgsMapToPixel &mtp = context.mapToPixel(); - QgsLineString pts; + QVector< double > pointsX; + QVector< double > pointsY; + QVector< double > pointsZ; // apply clipping for large lines to achieve a better rendering performance if ( clipToExtent && nPoints > 1 && !( context.flags() & Qgis::RenderContextFlag::ApplyClipAfterReprojection ) ) @@ -104,30 +106,51 @@ QPolygonF QgsSymbol::_getLineString3d( QgsRenderContext &context, const QgsCurve const double cw = e.width() / 10; const double ch = e.height() / 10; const QgsBox3d clipRect( e.xMinimum() - cw, e.yMinimum() - ch, -HUGE_VAL, e.xMaximum() + cw, e.yMaximum() + ch, HUGE_VAL ); // TODO also need to be clipped according to z axis - pts = QgsClipper::clipped3dLine( curve, clipRect ); + + const QgsLineString *lineString = nullptr; + if ( const QgsLineString *ls = qgsgeometry_cast< const QgsLineString * >( &curve ) ) + { + lineString = ls; + } + else + { + std::unique_ptr< QgsLineString > segmentized; + segmentized.reset( qgsgeometry_cast< QgsLineString * >( curve.segmentize( ) ) ); + lineString = segmentized.get(); + } + + QgsClipper::clipped3dLine( lineString->xVector(), lineString->yVector(), lineString->zVector(), pointsX, pointsY, pointsZ, clipRect ); } else { // clone... - const QgsLineString *tmpLine = reinterpret_cast( &curve ); - if ( dynamic_cast( &curve ) ) + if ( const QgsLineString *ls = qgsgeometry_cast( &curve ) ) { - pts.setPoints( curve.numPoints(), tmpLine->xData(), tmpLine->yData(), tmpLine->zData(), tmpLine->mData() ); + pointsX = ls->xVector(); + pointsY = ls->yVector(); + pointsZ = ls->zVector(); } else { - QgsPointSequence seq; - curve.points( seq ); - pts.setPoints( seq ); + std::unique_ptr< QgsLineString > segmentized; + segmentized.reset( qgsgeometry_cast< QgsLineString * >( curve.segmentize( ) ) ); + + pointsX = segmentized->xVector(); + pointsY = segmentized->yVector(); + pointsZ = segmentized->zVector(); } } - // transform the QPolygonF to screen coordinates + // transform the points to screen coordinates if ( ct.isValid() ) { + //create x, y arrays + const int nVertices = pointsX.size(); + + QString err; try { - pts.transform( ct, Qgis::TransformDirection::Forward, true ); + ct.transformCoords( nVertices, pointsX.data(), pointsY.data(), pointsZ.data(), Qgis::TransformDirection::Forward ); } catch ( QgsCsException & ) { @@ -136,17 +159,36 @@ QPolygonF QgsSymbol::_getLineString3d( QgsRenderContext &context, const QgsCurve } // remove non-finite points, e.g. infinite or NaN points caused by reprojecting errors - pts.filterVertices( []( const QgsPoint & point ) { - if ( point.is3D() ) + const int size = pointsX.size(); + + const double *xIn = pointsX.data(); + const double *yIn = pointsY.data(); + const double *zIn = pointsZ.data(); + double *xOut = pointsX.data(); + double *yOut = pointsY.data(); + double *zOut = pointsZ.data(); + int outSize = 0; + for ( int i = 0; i < size; ++i ) { - return std::isfinite( point.x() ) && std::isfinite( point.y() ) && std::isfinite( point.z() ); + if ( std::isfinite( *xIn ) && std::isfinite( *yIn ) && std::isfinite( *zIn ) ) + { + *xOut++ = *xIn++; + *yOut++ = *yIn++; + *zOut++ = *zIn++; + outSize++; + } + else + { + xIn++; + yIn++; + zIn++; + } } - else - { - return std::isfinite( point.x() ) && std::isfinite( point.y() ); - } - } ); + pointsX.resize( outSize ); + pointsY.resize( outSize ); + pointsZ.resize( outSize ); + } if ( clipToExtent && nPoints > 1 && context.flags() & Qgis::RenderContextFlag::ApplyClipAfterReprojection ) { @@ -155,14 +197,27 @@ QPolygonF QgsSymbol::_getLineString3d( QgsRenderContext &context, const QgsCurve const double cw = e.width() / 10; const double ch = e.height() / 10; const QgsBox3d clipRect( e.xMinimum() - cw, e.yMinimum() - ch, -HUGE_VAL, e.xMaximum() + cw, e.yMaximum() + ch, HUGE_VAL ); // TODO also need to be clipped according to z axis - pts = QgsClipper::clipped3dLine( pts, clipRect ); + + QVector< double > tempX; + QVector< double > tempY; + QVector< double > tempZ; + QgsClipper::clipped3dLine( pointsX, pointsY, pointsZ, tempX, tempY, tempZ, clipRect ); + pointsX = tempX; + pointsY = tempY; + pointsZ = tempZ; } - QPolygonF out = pts.asQPolygonF(); - QPointF *ptr = out.data(); - for ( int i = 0; i < out.size(); ++i, ++ptr ) + const int polygonSize = pointsX.size(); + QPolygonF out( polygonSize ); + const double *x = pointsX.constData(); + const double *y = pointsY.constData(); + QPointF *dest = out.data(); + for ( int i = 0; i < polygonSize; ++i ) { - mtp.transformInPlace( ptr->rx(), ptr->ry() ); + double screenX = *x++; + double screenY = *y++; + mtp.transformInPlace( screenX, screenY ); + *dest++ = QPointF( screenX, screenY ); } return out; @@ -243,40 +298,79 @@ QPolygonF QgsSymbol::_getPolygonRing3d( QgsRenderContext &context, const QgsCurv const QgsCoordinateTransform ct = context.coordinateTransform(); const QgsMapToPixel &mtp = context.mapToPixel(); - // clone... - QgsPointSequence seq; - curve.points( seq ); - QgsLineString poly( seq ); + QVector< double > pointsX; + QVector< double > pointsY; + QVector< double > pointsZ; if ( curve.numPoints() < 1 ) return QPolygonF(); + bool reverseRing = false; if ( correctRingOrientation ) { // ensure consistent polygon ring orientation if ( ( isExteriorRing && curve.orientation() != Qgis::AngularDirection::Clockwise ) || ( !isExteriorRing && curve.orientation() != Qgis::AngularDirection::CounterClockwise ) ) { - poly.reversed()->points( seq ); - poly.setPoints( seq ); + reverseRing = true; } } //clip close to view extent, if needed - if ( clipToExtent && !( context.flags() & Qgis::RenderContextFlag::ApplyClipAfterReprojection ) && !context.extent().contains( poly.boundingBox() ) ) + if ( clipToExtent && !( context.flags() & Qgis::RenderContextFlag::ApplyClipAfterReprojection ) && !context.extent().contains( curve.boundingBox() ) ) { const QgsRectangle e = context.extent(); const double cw = e.width() / 10; const double ch = e.height() / 10; const QgsBox3d clipRect( e.xMinimum() - cw, e.yMinimum() - ch, -HUGE_VAL, e.xMaximum() + cw, e.yMaximum() + ch, HUGE_VAL ); // TODO also need to be clipped according to z axis - QgsClipper::trimPolygon( poly, clipRect ); + + const QgsLineString *lineString = nullptr; + if ( const QgsLineString *ls = qgsgeometry_cast< const QgsLineString * >( &curve ) ) + { + lineString = ls; + } + else + { + std::unique_ptr< QgsLineString > segmentized; + segmentized.reset( qgsgeometry_cast< QgsLineString * >( curve.segmentize( ) ) ); + lineString = segmentized.get(); + } + + QgsClipper::clipped3dLine( lineString->xVector(), lineString->yVector(), lineString->zVector(), pointsX, pointsY, pointsZ, clipRect ); + } + else + { + // clone... + if ( const QgsLineString *ls = qgsgeometry_cast( &curve ) ) + { + pointsX = ls->xVector(); + pointsY = ls->yVector(); + pointsZ = ls->zVector(); + } + else + { + std::unique_ptr< QgsLineString > segmentized; + segmentized.reset( qgsgeometry_cast< QgsLineString * >( curve.segmentize( ) ) ); + + pointsX = segmentized->xVector(); + pointsY = segmentized->yVector(); + pointsZ = segmentized->zVector(); + } + } + + if ( reverseRing ) + { + std::reverse( pointsX.begin(), pointsX.end() ); + std::reverse( pointsY.begin(), pointsY.end() ); + std::reverse( pointsZ.begin(), pointsZ.end() ); } //transform the QPolygonF to screen coordinates if ( ct.isValid() ) { + const int nVertices = pointsX.size(); try { - poly.transform( ct, Qgis::TransformDirection::Forward, true ); + ct.transformCoords( nVertices, pointsX.data(), pointsY.data(), pointsZ.data(), Qgis::TransformDirection::Forward ); } catch ( QgsCsException & ) { @@ -285,37 +379,68 @@ QPolygonF QgsSymbol::_getPolygonRing3d( QgsRenderContext &context, const QgsCurv } // remove non-finite points, e.g. infinite or NaN points caused by reprojecting errors - poly.filterVertices( []( const QgsPoint point ) { - if ( point.is3D() ) - { - return std::isfinite( point.x() ) && std::isfinite( point.y() ) && std::isfinite( point.z() ); - } - else - { - return std::isfinite( point.x() ) && std::isfinite( point.y() ); - } - } ); + const int size = pointsX.size(); + const double *xIn = pointsX.data(); + const double *yIn = pointsY.data(); + const double *zIn = pointsZ.data(); + double *xOut = pointsX.data(); + double *yOut = pointsY.data(); + double *zOut = pointsZ.data(); + int outSize = 0; + for ( int i = 0; i < size; ++i ) + { + if ( std::isfinite( *xIn ) && std::isfinite( *yIn ) && std::isfinite( *zIn ) ) + { + *xOut++ = *xIn++; + *yOut++ = *yIn++; + *zOut++ = *zIn++; + outSize++; + } + else + { + xIn++; + yIn++; + zIn++; + } + } + pointsX.resize( outSize ); + pointsY.resize( outSize ); + pointsZ.resize( outSize ); + } - if ( clipToExtent && context.flags() & Qgis::RenderContextFlag::ApplyClipAfterReprojection && !context.mapExtent().contains( poly.boundingBox() ) ) + if ( clipToExtent && context.flags() & Qgis::RenderContextFlag::ApplyClipAfterReprojection && !context.mapExtent().contains( curve.boundingBox() ) ) { // early clipping was not possible, so we have to apply it here after transformation const QgsRectangle e = context.mapExtent(); const double cw = e.width() / 10; const double ch = e.height() / 10; const QgsBox3d clipRect( e.xMinimum() - cw, e.yMinimum() - ch, -HUGE_VAL, e.xMaximum() + cw, e.yMaximum() + ch, HUGE_VAL ); // TODO also need to be clipped according to z axis - QgsClipper::trimPolygon( poly, clipRect ); + + QVector< double > tempX; + QVector< double > tempY; + QVector< double > tempZ; + QgsClipper::clipped3dLine( pointsX, pointsY, pointsZ, tempX, tempY, tempZ, clipRect ); + pointsX = tempX; + pointsY = tempY; + pointsZ = tempZ; } - QPolygonF out = poly.asQPolygonF(); - QPointF *ptr = out.data(); - for ( int i = 0; i < out.size(); ++i, ++ptr ) + const int polygonSize = pointsX.size(); + QPolygonF out( polygonSize ); + const double *x = pointsX.constData(); + const double *y = pointsY.constData(); + QPointF *dest = out.data(); + for ( int i = 0; i < polygonSize; ++i ) { - mtp.transformInPlace( ptr->rx(), ptr->ry() ); + double screenX = *x++; + double screenY = *y++; + mtp.transformInPlace( screenX, screenY ); + *dest++ = QPointF( screenX, screenY ); } - if ( !out.empty() && !poly.isClosed() ) + if ( !out.empty() && !out.isClosed() ) out << out.at( 0 ); return out; diff --git a/tests/src/core/testqgsclipper.cpp b/tests/src/core/testqgsclipper.cpp index 539c0286155..b110d8cd4d8 100644 --- a/tests/src/core/testqgsclipper.cpp +++ b/tests/src/core/testqgsclipper.cpp @@ -50,77 +50,77 @@ void TestQgsClipper::initTestCase() void TestQgsClipper::basicWithZ() { // QgsClipper is static only - - QgsLineString polygon; - polygon.addVertex( QgsPoint( 10.4, 20.5, 10.0 ) ); - polygon.addVertex( QgsPoint( 20.2, 30.2, 20.0 ) ); - QgsBox3d clipRect( 10, 10, 11, 25, 30, 19 ); - QgsClipper::trimPolygon( polygon, clipRect ); + QVector< double > x = { 10.4, 20.2 }; + QVector< double > y = {20.5, 30.2 }; + QVector< double > z = {10.0, 20.0 }; + QgsClipper::trimPolygon( x, y, z, clipRect ); // Check nothing sticks out. - QVERIFY( checkBoundingBox( polygon, clipRect ) ); + QVERIFY( checkBoundingBox( QgsLineString( x, y, z ), clipRect ) ); // Check that it didn't clip too much QgsBox3d clipRectInner( clipRect ); clipRectInner.scale( 0.999 ); - QVERIFY( ! checkBoundingBox( polygon, clipRectInner ) ); + QVERIFY( ! checkBoundingBox( QgsLineString( x, y, z ), clipRectInner ) ); // A more complex example - polygon.clear(); - polygon.addVertex( QgsPoint( 1.0, 9.0, 1.0 ) ); - polygon.addVertex( QgsPoint( 11.0, 11.0, 11.0 ) ); - polygon.addVertex( QgsPoint( 9.0, 1.0, 9.0 ) ); + x = { 1.0, 11.0, 9.0 }; + y = { 9.0, 11.0, 1.0 }; + z = { 1.0, 11.0, 9.0 }; clipRect = QgsBox3d( 0.0, 0.0, 0.0, 10.0, 10.0, 10.0 ); - QgsClipper::trimPolygon( polygon, clipRect ); + QgsClipper::trimPolygon( x, y, z, clipRect ); // We should have 5 vertices now? - QCOMPARE( polygon.numPoints(), 5 ); + QCOMPARE( x.size(), 5 ); + QCOMPARE( y.size(), 5 ); + QCOMPARE( z.size(), 5 ); // Check nothing sticks out. - QVERIFY( checkBoundingBox( polygon, clipRect ) ); + QVERIFY( checkBoundingBox( QgsLineString( x, y, z ), clipRect ) ); // Check that it didn't clip too much clipRectInner = clipRect; clipRectInner.scale( 0.999 ); - QVERIFY( ! checkBoundingBox( polygon, clipRectInner ) ); + QVERIFY( ! checkBoundingBox( QgsLineString( x, y, z ), clipRectInner ) ); } void TestQgsClipper::basicWithZInf() { // QgsClipper is static only - QgsLineString polygon; - polygon.addVertex( QgsPoint( 10.4, 20.5, 10.0 ) ); - polygon.addVertex( QgsPoint( 20.2, 30.2, 20.0 ) ); + QVector< double > x { 10.4, 20.2 }; + QVector< double > y { 20.5, 30.2 }; + QVector< double > z { 10.0, 20.0 }; QgsBox3d clipRect( 10, 10, -HUGE_VAL, 25, 30, HUGE_VAL ); - QgsClipper::trimPolygon( polygon, clipRect ); + QgsClipper::trimPolygon( x, y, z, clipRect ); // Check nothing sticks out. - QVERIFY( checkBoundingBox( polygon, clipRect ) ); + QVERIFY( checkBoundingBox( QgsLineString( x, y, z ), clipRect ) ); // Check that it didn't clip too much QgsBox3d clipRectInner( clipRect ); clipRectInner.scale( 0.999 ); - QVERIFY( ! checkBoundingBox( polygon, clipRectInner ) ); + QVERIFY( ! checkBoundingBox( QgsLineString( x, y, z ), clipRectInner ) ); // A more complex example - polygon.clear(); - polygon.addVertex( QgsPoint( 1.0, 9.0, 1.0 ) ); - polygon.addVertex( QgsPoint( 11.0, 11.0, 11.0 ) ); - polygon.addVertex( QgsPoint( 9.0, 1.0, 9.0 ) ); + x = { 1.0, 11.0, 9.0 }; + y = { 9.0, 11.0, 1.0 }; + z = { 1.0, 11.0, 9.0 }; clipRect = QgsBox3d( 0.0, 0.0, 0.0, 10.0, 10.0, 10.0 ); - QgsClipper::trimPolygon( polygon, clipRect ); + QgsClipper::trimPolygon( x, y, z, clipRect ); // We should have 5 vertices now? - QCOMPARE( polygon.numPoints(), 5 ); + QCOMPARE( x.size(), 5 ); + QCOMPARE( y.size(), 5 ); + QCOMPARE( z.size(), 5 ); // Check nothing sticks out. - QVERIFY( checkBoundingBox( polygon, clipRect ) ); + QVERIFY( checkBoundingBox( QgsLineString( x, y, z ), clipRect ) ); // Check that it didn't clip too much clipRectInner = clipRect; clipRectInner.scale( 0.999 ); - QVERIFY( ! checkBoundingBox( polygon, clipRectInner ) ); + QVERIFY( ! checkBoundingBox( QgsLineString( x, y, z ), clipRectInner ) ); } void TestQgsClipper::basic()