Skip to content
Permalink
Browse files

Make pole of inaccessibility calculation handle multipolygons

  • Loading branch information
nyalldawson committed Jul 16, 2017
1 parent 6487fbb commit 38a13ff5af23e8c56df656a880b4ccf4a584362d
Showing with 68 additions and 21 deletions.
  1. +61 −21 src/core/geometry/qgsinternalgeometryengine.cpp
  2. +7 −0 tests/src/core/testqgsgeometry.cpp
@@ -155,25 +155,10 @@ Cell *getCentroidCell( const QgsPolygonV2 *polygon )
return new Cell( x / area, y / area, 0.0, polygon );
}

///@endcond

QgsGeometry QgsInternalGeometryEngine::poleOfInaccessibility( double precision, double *distanceFromBoundary ) const
QgsPoint surfacePoleOfInaccessibility( const QgsSurface *surface, double precision, double &distanceFromBoundary )
{
if ( distanceFromBoundary )
*distanceFromBoundary = DBL_MAX;

if ( !mGeometry || mGeometry->isEmpty() )
return QgsGeometry();

if ( precision <= 0 )
return QgsGeometry();

const QgsSurface *surface = dynamic_cast< const QgsSurface * >( mGeometry );
if ( !surface )
return QgsGeometry();

std::unique_ptr< QgsPolygonV2 > segmentizedPoly;
const QgsPolygonV2 *polygon = dynamic_cast< const QgsPolygonV2 * >( mGeometry );
const QgsPolygonV2 *polygon = dynamic_cast< const QgsPolygonV2 * >( surface );
if ( !polygon )
{
segmentizedPoly.reset( static_cast< QgsPolygonV2 *>( surface->segmentize() ) );
@@ -187,7 +172,7 @@ QgsGeometry QgsInternalGeometryEngine::poleOfInaccessibility( double precision,
double cellSize = qMin( bounds.width(), bounds.height() );

if ( qgsDoubleNear( cellSize, 0.0 ) )
return QgsGeometry( new QgsPoint( bounds.xMinimum(), bounds.yMinimum() ) );
return QgsPoint( bounds.xMinimum(), bounds.yMinimum() );

double h = cellSize / 2.0;
std::priority_queue< Cell *, std::vector<Cell *>, GreaterThanByMax > cellQueue;
@@ -238,15 +223,70 @@ QgsGeometry QgsInternalGeometryEngine::poleOfInaccessibility( double precision,
cellQueue.push( new Cell( currentCell->x + h, currentCell->y + h, h, polygon ) );
}

distanceFromBoundary = bestCell->d;

return QgsPoint( bestCell->x, bestCell->y );
}

///@endcond

QgsGeometry QgsInternalGeometryEngine::poleOfInaccessibility( double precision, double *distanceFromBoundary ) const
{
if ( distanceFromBoundary )
*distanceFromBoundary = bestCell->d;
*distanceFromBoundary = DBL_MAX;

return QgsGeometry( new QgsPoint( bestCell->x, bestCell->y ) );
if ( !mGeometry || mGeometry->isEmpty() )
return QgsGeometry();

if ( precision <= 0 )
return QgsGeometry();

if ( const QgsGeometryCollection *gc = dynamic_cast< const QgsGeometryCollection *>( mGeometry ) )
{
int numGeom = gc->numGeometries();
double maxDist = 0;
QgsPoint bestPoint;
bool found = false;
for ( int i = 0; i < numGeom; ++i )
{
const QgsSurface *surface = dynamic_cast< const QgsSurface * >( gc->geometryN( i ) );
if ( !surface )
continue;

found = true;
double dist = DBL_MAX;
QgsPoint p = surfacePoleOfInaccessibility( surface, precision, dist );
if ( dist > maxDist )
{
maxDist = dist;
bestPoint = p;
}
}

if ( !found )
return QgsGeometry();

if ( distanceFromBoundary )
*distanceFromBoundary = maxDist;
return QgsGeometry( new QgsPoint( bestPoint ) );
}
else
{
const QgsSurface *surface = dynamic_cast< const QgsSurface * >( mGeometry );
if ( !surface )
return QgsGeometry();

double dist = DBL_MAX;
QgsPoint p = surfacePoleOfInaccessibility( surface, precision, dist );
if ( distanceFromBoundary )
*distanceFromBoundary = dist;
return QgsGeometry( new QgsPoint( p ) );
}
}


// helpers for orthogonalize
// adapted from original code in potlach/id osm editor
// adapted from original code in potlatch/id osm editor

bool dotProductWithinAngleTolerance( double dotProduct, double lowerThreshold, double upperThreshold )
{
@@ -5451,6 +5451,13 @@ void TestQgsGeometry::poleOfInaccessibility()
point = curved.poleOfInaccessibility( 0.01 ).asPoint();
QGSCOMPARENEAR( point.x(), -0.4324, 0.0001 );
QGSCOMPARENEAR( point.y(), -0.2434, 0.0001 );

// multipolygon
QgsGeometry multiPoly = QgsGeometry::fromWkt( QStringLiteral( "MultiPolygon (((0 0, 10 0, 10 10, 0 10, 0 0)),((30 30, 50 30, 50 60, 30 60, 30 30)))" ) );
point = multiPoly.poleOfInaccessibility( 0.01, &distance ).asPoint();
QGSCOMPARENEAR( point.x(), 40, 0.0001 );
QGSCOMPARENEAR( point.y(), 45, 0.0001 );
QGSCOMPARENEAR( distance, 10.0, 0.00001 );
}

void TestQgsGeometry::makeValid()

0 comments on commit 38a13ff

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