Skip to content

Commit

Permalink
Browse files Browse the repository at this point in the history
Add QgsGeometry::isAxisParallelRectangle for checking whether
a geometry is a axis-parallel rectangle

Credit to @stefanuhrig
  • Loading branch information
nyalldawson committed Apr 23, 2021
1 parent 6470771 commit 9c4375b
Show file tree
Hide file tree
Showing 6 changed files with 226 additions and 0 deletions.
15 changes: 15 additions & 0 deletions python/core/auto_generated/geometry/qgsgeometry.sip.in
Expand Up @@ -353,6 +353,21 @@ Uses GEOS library for the test.
for checking anomalies in polygon geometries one can use :py:func:`~QgsGeometry.isGeosValid`.

.. versionadded:: 3.0
%End

bool isAxisParallelRectangle( double maximumDeviation, bool simpleRectanglesOnly = false ) const;
%Docstring
Returns ``True`` if the geometry is a polygon that is almost an axis-parallel rectangle.

The ``maximumDeviation`` argument specifes the maximum angle (in degrees) that the polygon edges
are allowed to deviate from axis parallel lines.

By default the check will permit polygons with more than 4 edges, so long as the overall shape of
the polygon is an axis-parallel rectangle (i.e. it is tolerant to rectangles with additional vertices
added along the rectangle sides). If ``simpleRectanglesOnly`` is set to ``True`` then the method will
only return ``True`` if the geometry is a simple rectangle consisting of 4 edges.

.. versionadded:: 3.20
%End

double area() const;
Expand Down
9 changes: 9 additions & 0 deletions src/core/geometry/qgsgeometry.cpp
Expand Up @@ -2737,6 +2737,15 @@ bool QgsGeometry::isSimple() const
return geos.isSimple( &mLastError );
}

bool QgsGeometry::isAxisParallelRectangle( double maximumDeviation, bool simpleRectanglesOnly ) const
{
if ( !d->geometry )
return false;

QgsInternalGeometryEngine engine( *this );
return engine.isAxisParallelRectangle( maximumDeviation, simpleRectanglesOnly );
}

bool QgsGeometry::isGeosEqual( const QgsGeometry &g ) const
{
if ( !d->geometry || !g.d->geometry )
Expand Down
15 changes: 15 additions & 0 deletions src/core/geometry/qgsgeometry.h
Expand Up @@ -389,6 +389,21 @@ class CORE_EXPORT QgsGeometry
*/
bool isSimple() const;

/**
* Returns TRUE if the geometry is a polygon that is almost an axis-parallel rectangle.
*
* The \a maximumDeviation argument specifes the maximum angle (in degrees) that the polygon edges
* are allowed to deviate from axis parallel lines.
*
* By default the check will permit polygons with more than 4 edges, so long as the overall shape of
* the polygon is an axis-parallel rectangle (i.e. it is tolerant to rectangles with additional vertices
* added along the rectangle sides). If \a simpleRectanglesOnly is set to TRUE then the method will
* only return TRUE if the geometry is a simple rectangle consisting of 4 edges.
*
* \since QGIS 3.20
*/
bool isAxisParallelRectangle( double maximumDeviation, bool simpleRectanglesOnly = false ) const;

/**
* Returns the planar, 2-dimensional area of the geometry.
*
Expand Down
133 changes: 133 additions & 0 deletions src/core/geometry/qgsinternalgeometryengine.cpp
Expand Up @@ -46,6 +46,139 @@ QString QgsInternalGeometryEngine::lastError() const
return mLastError;
}


enum class Direction
{
Up,
Right,
Down,
Left,
None
};

/**
* Determines the direction of an edge from p1 to p2. maxDev is the tangent of
* the maximum allowed edge deviation angle. If the edge deviates more than
* the allowed angle, Direction::None will be returned.
*/
Direction getEdgeDirection( const QgsPoint &p1, const QgsPoint &p2, double maxDev )
{
double dx = p2.x() - p1.x();
double dy = p2.y() - p1.y();
if ( ( dx == 0.0 ) && ( dy == 0.0 ) )
return Direction::None;
if ( fabs( dx ) >= fabs( dy ) )
{
double dev = fabs( dy ) / fabs( dx );
if ( dev > maxDev )
return Direction::None;
return dx > 0.0 ? Direction::Right : Direction::Left;
}
else
{
double dev = fabs( dx ) / fabs( dy );
if ( dev > maxDev )
return Direction::None;
return dy > 0.0 ? Direction::Up : Direction::Down;
}
}

/**
* Checks whether the polygon consists of four nearly axis-parallel sides. All
* consecutive edges having the same direction are considered to belong to the
* same side.
*/
std::pair<bool, std::array<Direction, 4>> getEdgeDirections( const QgsPolygon *g, double maxDev )
{
std::pair<bool, std::array<Direction, 4>> ret = { false, { Direction::None, Direction::None, Direction::None, Direction::None } };
// The polygon might start in the middle of a side. Hence, we need a fifth
// direction to record the beginning of the side when we went around the
// polygon.
std::array<Direction, 5> dirs;

int idx = 0;
QgsAbstractGeometry::vertex_iterator previous = g->vertices_begin();
QgsAbstractGeometry::vertex_iterator current = previous;
++current;
QgsAbstractGeometry::vertex_iterator end = g->vertices_end();
while ( current != end )
{
Direction dir = getEdgeDirection( *previous, *current, maxDev );
if ( dir == Direction::None )
return ret;
if ( idx == 0 )
{
dirs[0] = dir;
++idx;
}
else if ( dir != dirs[idx - 1] )
{
if ( idx == 5 )
return ret;
dirs[idx] = dir;
++idx;
}
previous = current;
++current;
}
ret.first = ( idx == 5 ) ? ( dirs[0] == dirs[4] ) : ( idx == 4 );
std::copy( dirs.begin(), dirs.begin() + 4, ret.second.begin() );
return ret;
}

bool matchesOrientation( std::array<Direction, 4> dirs, std::array<Direction, 4> oriented )
{
int idx = std::find( oriented.begin(), oriented.end(), dirs[0] ) - oriented.begin();
for ( int i = 1; i < 4; ++i )
{
if ( dirs[i] != oriented[( idx + i ) % 4] )
return false;
}
return true;
}

/**
* Checks whether the 4 directions in dirs make up a clockwise rectangle.
*/
bool isClockwise( std::array<Direction, 4> dirs )
{
const std::array<Direction, 4> cwdirs = { Direction::Up, Direction::Right, Direction::Down, Direction::Left };
return matchesOrientation( dirs, cwdirs );
}

/**
* Checks whether the 4 directions in dirs make up a counter-clockwise
* rectangle.
*/
bool isCounterClockwise( std::array<Direction, 4> dirs )
{
const std::array<Direction, 4> ccwdirs = { Direction::Right, Direction::Up, Direction::Left, Direction::Down };
return matchesOrientation( dirs, ccwdirs );
}


bool QgsInternalGeometryEngine::isAxisParallelRectangle( double maximumDeviation, bool simpleRectanglesOnly ) const
{
if ( QgsWkbTypes::flatType( mGeometry->wkbType() ) != QgsWkbTypes::Polygon )
return false;

const QgsPolygon *polygon = qgsgeometry_cast< const QgsPolygon * >( mGeometry );
if ( !polygon->exteriorRing() || polygon->numInteriorRings() > 0 )
return false;

const int vertexCount = polygon->exteriorRing()->numPoints();
if ( vertexCount < 4 )
return false;
else if ( simpleRectanglesOnly && ( vertexCount != 5 || !polygon->exteriorRing()->isClosed() ) )
return false;

bool found4Dirs;
std::array<Direction, 4> dirs;
std::tie( found4Dirs, dirs ) = getEdgeDirections( polygon, std::tan( maximumDeviation * M_PI / 180 ) );

return found4Dirs && ( isCounterClockwise( dirs ) || isClockwise( dirs ) );
}

/***************************************************************************
* This class is considered CRITICAL and any change MUST be accompanied with
* full unit tests.
Expand Down
15 changes: 15 additions & 0 deletions src/core/geometry/qgsinternalgeometryengine.h
Expand Up @@ -56,6 +56,21 @@ class QgsInternalGeometryEngine
*/
QString lastError() const;

/**
* Returns TRUE if the geometry is a polygon that is almost an axis-parallel rectangle.
*
* The \a maximumDeviation argument specifes the maximum angle (in degrees) that the polygon edges
* are allowed to deviate from axis parallel lines.
*
* By default the check will permit polygons with more than 4 edges, so long as the overall shape of
* the polygon is an axis-parallel rectangle (i.e. it is tolerant to rectangles with additional vertices
* added along the rectangle sides). If \a simpleRectanglesOnly is set to TRUE then the method will
* only return TRUE if the geometry is a simple rectangle consisting of 4 edges.
*
* \since QGIS 3.20
*/
bool isAxisParallelRectangle( double maximumDeviation, bool simpleRectanglesOnly = false ) const;

/**
* Will extrude a line or (segmentized) curve by a given offset and return a polygon
* representation of it.
Expand Down
39 changes: 39 additions & 0 deletions tests/src/python/test_qgsgeometry.py
Expand Up @@ -6011,6 +6011,45 @@ def testGeosCrash(self):
# test we don't crash when geos returns a point geometry with no points
QgsGeometry.fromWkt('Polygon ((0 0, 1 1, 1 0, 0 0))').intersection(QgsGeometry.fromWkt('Point (42 0)')).isNull()

def testIsRectangle(self):
"""
Test checking if geometries are rectangles
"""
# non polygons
self.assertFalse(QgsGeometry().isAxisParallelRectangle(0))
self.assertFalse(QgsGeometry.fromWkt('Point(0 1)').isAxisParallelRectangle(0))
self.assertFalse(QgsGeometry.fromWkt('LineString(0 1, 1 2)').isAxisParallelRectangle(0))

self.assertFalse(QgsGeometry.fromWkt('Polygon((0 0, 0 1, 1 1, 0 0))').isAxisParallelRectangle(0))
self.assertFalse(QgsGeometry.fromWkt('Polygon((0 0, 0 1, 1 1, 0 0))').isAxisParallelRectangle(0, True))
self.assertFalse(QgsGeometry.fromWkt('Polygon((0 0, 0 1, 1 1, 0 1))').isAxisParallelRectangle(0))
self.assertFalse(QgsGeometry.fromWkt('Polygon((0 0, 0 1, 1 1, 0 1))').isAxisParallelRectangle(0, True))
self.assertFalse(QgsGeometry.fromWkt('Polygon(())').isAxisParallelRectangle(0))

self.assertTrue(QgsGeometry.fromWkt('Polygon((0 0, 1 0, 1 1, 0 1, 0 0))').isAxisParallelRectangle(0))
self.assertTrue(QgsGeometry.fromWkt('Polygon((0 0, 1 0, 1 1, 0 1, 0 0))').isAxisParallelRectangle(0, True))
self.assertTrue(QgsGeometry.fromWkt('Polygon((0 0, 0 1, 1 1, 1 0, 0 0))').isAxisParallelRectangle(0))
self.assertTrue(QgsGeometry.fromWkt('Polygon((0 0, 0 1, 1 1, 1 0, 0 0))').isAxisParallelRectangle(0, True))
# with rings
self.assertFalse(QgsGeometry.fromWkt('Polygon((0 0, 1 0, 1 1, 0 1, 0 0), (0.1 0.1, 0.2 0.1, 0.2 0.2, 0.1 0.2, 0.1 0.1))').isAxisParallelRectangle(0))
# densified
self.assertTrue(QgsGeometry.fromWkt('Polygon((0 0, 0.5 0.0, 1 0, 1 1, 0 1, 0 0))').densifyByCount(5).isAxisParallelRectangle(0))
# not a simple rectangle
self.assertFalse(QgsGeometry.fromWkt('Polygon((0 0, 0.5 0.0, 1 0, 1 1, 0 1, 0 0))').densifyByCount(5).isAxisParallelRectangle(0, True))

# starting mid way through a side
self.assertTrue(QgsGeometry.fromWkt('Polygon((0.5 0, 1 0, 1 1, 0 1, 0 0, 0.5 0))').densifyByCount(
5).isAxisParallelRectangle(0))
self.assertFalse(QgsGeometry.fromWkt('Polygon((0.5 0, 1 0, 1 1, 0 1, 0 0, 0.5 0))').densifyByCount(
5).isAxisParallelRectangle(0, True))

# with tolerance
self.assertFalse(QgsGeometry.fromWkt('Polygon((0 0, 1 0.001, 1 1, 0 1, 0 0))').isAxisParallelRectangle(0))
self.assertTrue(QgsGeometry.fromWkt('Polygon((0 0, 1 0.001, 1 1, 0 1, 0 0))').isAxisParallelRectangle(1))
self.assertFalse(QgsGeometry.fromWkt('Polygon((0 0, 1 0.1, 1 1, 0 1, 0 0))').isAxisParallelRectangle(1))
self.assertTrue(QgsGeometry.fromWkt('Polygon((0 0, 1 0.1, 1 1, 0 1, 0 0))').isAxisParallelRectangle(10))
self.assertTrue(QgsGeometry.fromWkt('Polygon((0 0, 1 0.1, 1 1, 0 1, 0 0))').densifyByCount(5).isAxisParallelRectangle(10))

def testTransformWithClass(self):
"""
Test transforming using a python based class (most of the tests for this are in c++, this is just
Expand Down

0 comments on commit 9c4375b

Please sign in to comment.