sharpetronics/QGIS
1
0
Fork 0
mirror of https://github.com/qgis/QGIS.git synced 2025-10-06 00:07:29 -04:00
Code Issues Packages Projects Releases Wiki Activity

Optimise 3d clipping to work directly with coordinate arrays, so that

Browse Source 
we avoid the (signficant) overhead of point conversion and sequence
allocations
This commit is contained in:
Nyall Dawson 2022-05-25 15:20:35 +10:00
parent 54ab1dc4fe
commit f439a649d6
7 changed files with 430 additions and 203 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

5
python/core/auto_generated/geometry/qgslinestring.sip.in
View File

@ -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.

20
python/core/auto_generated/qgsclipper.sip.in
View File

@ -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

41
src/core/geometry/qgslinestring.h
View File

@ -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
/**

115
src/core/qgsclipper.cpp
View File

@ -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<double> &zIn, QVector<double> &x, QVector<double> &y, QVector<double> &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<double> &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;
}
}

169
src/core/qgsclipper.h
View File

@ -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<double> &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<double> &inX, const QVector<double> &inY, const QVector< double > &inZ,
QVector<double> &outX, QVector<double> &outY, QVector<double> &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<double> &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<double> &x, QVector<double> &y, QVector<double> &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<double> &inX, const QVector<double> &inY, const QVector< double > &inZ,
QVector<double> &outX, QVector<double> &outY, QVector<double> &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 )

223
src/core/symbology/qgssymbol.cpp
View File

@ -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<const QgsLineString *>( &curve );
if ( dynamic_cast<const QgsLineString *>( &curve ) )
if ( const QgsLineString *ls = qgsgeometry_cast<const QgsLineString *>( &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<const QgsLineString *>( &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;

60
tests/src/core/testqgsclipper.cpp
View File

@ -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()