Skip to content

Commit

Permalink
Optimize centroid calculation rubber sheet (#5627)
Browse files Browse the repository at this point in the history 
* Calculate the centroid of an element as the average of all points in the element instead of converting to GEOS::GEOM objects
* Integrate ElementGeometryUtils::calculateElementCentroid into CentroidDistanceExtractor
* Optimized ElementGeometryUtils::calculateLength
* Copyright
* Sonar findings
  • Loading branch information
@bmarchant
bmarchant authored Apr 19, 2023
1 parent 11d262c commit d4276ba
Show file tree
Hide file tree
Showing 8 changed files with 604 additions and 388 deletions.
53 changes: 52 additions & 1 deletion hoot-core-test/src/test/cpp/hoot/core/elements/ElementGeometryUtilsTest.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -22,7 +22,7 @@
* This will properly maintain the copyright information. Maxar
* copyrights will be updated automatically.
*
* @copyright Copyright (C) 2020, 2021, 2022 Maxar (http://www.maxar.com/)
* @copyright Copyright (C) 2020-2023 Maxar (http://www.maxar.com/)
*/

// Hoot
Expand All @@ -43,6 +43,7 @@ class ElementGeometryUtilsTest : public HootTestFixture
CPPUNIT_TEST(haveRelationshipTest);
CPPUNIT_TEST(haveRelationshipNullInputTest);
CPPUNIT_TEST(geometricRelationshipToStringTest);
CPPUNIT_TEST(calculateElementCentroid);
CPPUNIT_TEST_SUITE_END();

public:
Expand Down Expand Up @@ -162,6 +163,56 @@ class ElementGeometryUtilsTest : public HootTestFixture
}
CPPUNIT_ASSERT_EQUAL(QString("Unknown type.").toStdString(), exceptionMsg.toStdString());
}

void calculateElementCentroid()
{
// Setup nodes
std::vector<NodePtr> nodes(
{
std::make_shared<Node>(Status::Unknown1, 1, 9.0, 9.0),
std::make_shared<Node>(Status::Unknown1, 2, 11.0, 9.0),
std::make_shared<Node>(Status::Unknown1, 3, 11.0, 11.0),
std::make_shared<Node>(Status::Unknown1, 4, 9.0, 11.0),
std::make_shared<Node>(Status::Unknown1, 5, 8.0, 8.0),
std::make_shared<Node>(Status::Unknown1, 6, 16.0, 8.0),
std::make_shared<Node>(Status::Unknown1, 7, 16.0, 16.0),
std::make_shared<Node>(Status::Unknown1, 8, 8.0, 16.0),
});
// Setup ways
std::vector<WayPtr> ways(
{
std::make_shared<Way>(Status::Unknown1, 1),
std::make_shared<Way>(Status::Unknown1, 2),
std::make_shared<Way>(Status::Unknown1, 3),
std::make_shared<Way>(Status::Unknown1, 4)
});
ways[0]->addNodes({ 1, 2, 3, 4 });
ways[1]->addNodes({ 1, 2, 3, 4, 1 });
ways[2]->addNodes({ 5, 6, 7, 8 });
ways[3]->addNodes({ 5, 6, 7, 8, 5 });
// Setup relation
RelationPtr relation = std::make_shared<Relation>(Status::Unknown1, 1);
relation->addElement("inner", ElementType::Way, 2);
relation->addElement("outer", ElementType::Way, 4);
// Add all elements to the map
OsmMapPtr map = std::make_shared<OsmMap>();
map->addElements(nodes.begin(), nodes.end());
map->addElements(ways.begin(), ways.end());
map->addRelation(relation);
// Validate the centroids of each type
geos::geom::Coordinate n = ElementGeometryUtils::calculateElementCentroid(ElementId(ElementType::Node, 1), map);
CPPUNIT_ASSERT_EQUAL(geos::geom::Coordinate(9.0, 9.0), n);
geos::geom::Coordinate w1 = ElementGeometryUtils::calculateElementCentroid(ElementId(ElementType::Way, 1), map);
CPPUNIT_ASSERT_EQUAL(geos::geom::Coordinate(10.0, 10.0), w1);
geos::geom::Coordinate w2 = ElementGeometryUtils::calculateElementCentroid(ElementId(ElementType::Way, 2), map);
CPPUNIT_ASSERT_EQUAL(geos::geom::Coordinate(10.0, 10.0), w2);
geos::geom::Coordinate w3 = ElementGeometryUtils::calculateElementCentroid(ElementId(ElementType::Way, 3), map);
CPPUNIT_ASSERT_EQUAL(geos::geom::Coordinate(12.0, 12.0), w3);
geos::geom::Coordinate w4 = ElementGeometryUtils::calculateElementCentroid(ElementId(ElementType::Way, 4), map);
CPPUNIT_ASSERT_EQUAL(geos::geom::Coordinate(12.0, 12.0), w4);
geos::geom::Coordinate r = ElementGeometryUtils::calculateElementCentroid(ElementId(ElementType::Relation, 1), map);
CPPUNIT_ASSERT_EQUAL(geos::geom::Coordinate(11.0, 11.0), r);
}
};

CPPUNIT_TEST_SUITE_NAMED_REGISTRATION(ElementGeometryUtilsTest, "quick");
Expand Down
21 changes: 7 additions & 14 deletions hoot-core/src/main/cpp/hoot/core/algorithms/extractors/CentroidDistanceExtractor.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -23,11 +23,12 @@
* copyrights will be updated automatically.
*
* @copyright Copyright (C) 2005 VividSolutions (http://www.vividsolutions.com/)
* @copyright Copyright (C) 2015, 2017, 2018, 2019, 2020, 2021, 2022 Maxar (http://www.maxar.com/)
* @copyright Copyright (C) 2015-2023 Maxar (http://www.maxar.com/)
*/
#include "CentroidDistanceExtractor.h"

// hoot
#include <hoot/core/elements/ElementGeometryUtils.h>
#include <hoot/core/geometry/ElementToGeometryConverter.h>
#include <hoot/core/geometry/GeometryUtils.h>
#include <hoot/core/util/Factory.h>
Expand All @@ -43,22 +44,14 @@ HOOT_FACTORY_REGISTER(FeatureExtractor, CentroidDistanceExtractor)
double CentroidDistanceExtractor::distance(const OsmMap& map, const std::shared_ptr<const Element>& target,
const std::shared_ptr<const Element>& candidate) const
{
ElementToGeometryConverter ec(map.shared_from_this());
std::shared_ptr<Geometry> g1 = ec.convertToGeometry(target);
std::shared_ptr<Geometry> g2 = ec.convertToGeometry(candidate);
ConstOsmMapPtr pmap = map.shared_from_this();
Coordinate tc = ElementGeometryUtils::calculateElementCentroid(target->getElementId(), pmap);
Coordinate cc = ElementGeometryUtils::calculateElementCentroid(candidate->getElementId(), pmap);

if (g1->isEmpty() || g2->isEmpty())
return -1;

g1.reset(GeometryUtils::validateGeometry(g1.get()));
g2.reset(GeometryUtils::validateGeometry(g2.get()));
std::shared_ptr<Point> tc(g1->getCentroid());
std::shared_ptr<Point> cc(g2->getCentroid());

if (tc.get() == nullptr || cc.get() == nullptr)
if (tc == Coordinate::getNull() || cc == Coordinate::getNull())
return nullValue();

return tc->distance(cc.get());
return tc.distance(cc);
}

}
26 changes: 6 additions & 20 deletions hoot-core/src/main/cpp/hoot/core/algorithms/rubber-sheet/RubberSheet.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -36,6 +36,7 @@
#include <hoot/core/criterion/NotCriterion.h>
#include <hoot/core/criterion/PolygonCriterion.h>
#include <hoot/core/criterion/WayNodeCriterion.h>
#include <hoot/core/elements/ElementGeometryUtils.h>
#include <hoot/core/elements/MapProjector.h>
#include <hoot/core/elements/NodeToWayMap.h>
#include <hoot/core/geometry/ElementToGeometryConverter.h>
Expand Down Expand Up @@ -624,30 +625,15 @@ bool RubberSheet::_findBuildingTies()
continue;
usedIds.insert(eid1);

// Find the center point of each element as the tie-point
ElementPtr e1 = _map->getElement(eid1);
ElementPtr e2 = _map->getElement(eid2);
Coordinate c1 = ElementGeometryUtils::calculateElementCentroid(eid1, _map);
Coordinate c2 = ElementGeometryUtils::calculateElementCentroid(eid2, _map);

if (!e1 || !e2)
continue;

std::shared_ptr<Geometry> g1 = ec.convertToGeometry(e1);
std::shared_ptr<Geometry> g2 = ec.convertToGeometry(e2);

if (g1->isEmpty() || g2->isEmpty())
continue;

g1.reset(GeometryUtils::validateGeometry(g1.get()));
g2.reset(GeometryUtils::validateGeometry(g2.get()));
std::shared_ptr<Point> c1(g1->getCentroid());
std::shared_ptr<Point> c2(g2->getCentroid());

if (!c1 || !c2)
if (c1 == Coordinate::getNull() || c2 == Coordinate::getNull())
continue;

Tie t;
t.c1 = geos::geom::Coordinate(c1->getX(), c1->getY());
t.c2 = geos::geom::Coordinate(c2->getX(), c2->getY());
t.c1 = c1;
t.c2 = c2;
LOG_VART(t.c1);
LOG_VART(t.c2);
_ties.push_back(t);
Expand Down
194 changes: 182 additions & 12 deletions hoot-core/src/main/cpp/hoot/core/elements/ElementGeometryUtils.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -22,7 +22,7 @@
* This will properly maintain the copyright information. Maxar
* copyrights will be updated automatically.
*
* @copyright Copyright (C) 2020, 2021, 2022 Maxar (http://www.maxar.com/)
* @copyright Copyright (C) 2020-2023 Maxar (http://www.maxar.com/)
*/

#include "ElementGeometryUtils.h"
Expand All @@ -39,6 +39,8 @@
// GEOS
#include <geos/util/TopologyException.h>

using namespace geos::geom;

namespace hoot
{

Expand All @@ -50,8 +52,8 @@ bool ElementGeometryUtils::haveGeometricRelationship(const ConstElementPtr& elem
if (!element1 || !element2)
throw IllegalArgumentException("One of the input elements is null.");

std::shared_ptr<geos::geom::Geometry> geom1 = _getGeometry(element1, map);
std::shared_ptr<geos::geom::Geometry> geom2 = _getGeometry(element2, map);
std::shared_ptr<Geometry> geom1 = _getGeometry(element1, map);
std::shared_ptr<Geometry> geom2 = _getGeometry(element2, map);
bool haveRelationship = false;
if (geom1 && geom2)
{
Expand Down Expand Up @@ -97,12 +99,12 @@ bool ElementGeometryUtils::haveGeometricRelationship(const ConstElementPtr& elem
return haveRelationship;
}

std::shared_ptr<geos::geom::Geometry> ElementGeometryUtils::_getGeometry(const ConstElementPtr& element, ConstOsmMapPtr map)
std::shared_ptr<Geometry> ElementGeometryUtils::_getGeometry(const ConstElementPtr& element, ConstOsmMapPtr map)
{
if (!element)
throw IllegalArgumentException("The input element is null.");

std::shared_ptr<geos::geom::Geometry> newGeom;
std::shared_ptr<Geometry> newGeom;
QString errorMsg = "Feature passed to ElementGeometryUtils caused topology exception on conversion to a geometry: ";
try
{
Expand Down Expand Up @@ -147,15 +149,183 @@ Meters ElementGeometryUtils::calculateLength(const ConstElementPtr& e, const Con
// NOTE: Originally, we used isLinear. This was a bit too strict in that it wants evidence of
// being linear before the length is calculated. Conversely, this wants evidence that it is not
// linear before it will assume it doesn't have a length.
if (e->getElementType() != ElementType::Node && AreaCriterion().isSatisfied(e) == false)
Meters length = 0.0;
if (!AreaCriterion().isSatisfied(e))
{
switch (e->getElementType().getEnum())
{
case ElementType::Node:
// Nodes don't have a length
break;
case ElementType::Way:
{
ConstWayPtr w = std::static_pointer_cast<const Way>(e);
length = _calculateLength(w, constProvider);
break;
}
case ElementType::Relation:
{
ConstRelationPtr r = std::static_pointer_cast<const Relation>(e);
length = _calculateLength(r, constProvider);
break;
}
default:
break;
}
}
return length;
}

Coordinate ElementGeometryUtils::calculateElementCentroid(const ElementId& eid, const ConstOsmMapPtr& map)
{
if (!map)
return Coordinate::getNull();
if (!map->containsElement(eid))
return Coordinate::getNull();
double centroid_x = 0.0;
double centroid_y = 0.0;
double centroid_count = 0.0;
ConstElementPtr element = map->getElement(eid);
switch(eid.getType().getEnum())
{
case ElementType::Node:
{
ConstNodePtr node = std::dynamic_pointer_cast<const Node>(element);
if (!_getCentroidValues(node, map, centroid_x, centroid_y, centroid_count))
return Coordinate::getNull();
break;
}
case ElementType::Way:
{
ConstWayPtr way = std::dynamic_pointer_cast<const Way>(element);
if (!_getCentroidValues(way, map, centroid_x, centroid_y, centroid_count))
return Coordinate::getNull();
break;
}
case ElementType::Relation:
{
ConstRelationPtr r = std::dynamic_pointer_cast<const Relation>(element);
if (!_getCentroidValues(r, map, centroid_x, centroid_y, centroid_count))
return Coordinate::getNull();
break;
}
default:
return Coordinate::getNull();
}
if (centroid_count < 1.0)
return Coordinate::getNull();
// Divide by the node count to find centroid
centroid_x /= centroid_count;
centroid_y /= centroid_count;
return Coordinate(centroid_x, centroid_y);
}

bool ElementGeometryUtils::_getCentroidValues(const ConstNodePtr& node, const ConstOsmMapPtr& /*map*/, double& centroid_x, double& centroid_y, double& centroid_count)
{
if (!node)
return false;
centroid_x += node->getX();
centroid_y += node->getY();
centroid_count++;
return true;
}

bool ElementGeometryUtils::_getCentroidValues(const ConstWayPtr& way, const ConstOsmMapPtr& map, double& centroid_x, double& centroid_y, double& centroid_count)
{
if (!way || !map)
return false;
if (way->getNodeCount() < 1)
return false;
double x = 0.0;
double y = 0.0;
const std::vector<long> node_ids = way->getNodeIds();
size_t size = node_ids.size();
// Don't include the last ID if it is a closed way
if (node_ids[0] == node_ids[size - 1])
size--;
// Sum all way nodes
for (size_t index = 0; index < size; ++index)
{
long node_id = node_ids[index];
if (!map->containsNode(node_id))
return false;
ConstNodePtr node = map->getNode(node_id);
x += node->getX();
y += node->getY();
}
centroid_x += x;
centroid_y += y;
centroid_count += static_cast<double>(size);
return true;
}

bool ElementGeometryUtils::_getCentroidValues(const ConstRelationPtr& relation, const ConstOsmMapPtr& map, double& centroid_x, double& centroid_y, double& centroid_count)
{
if (!relation || !map)
return false;
double x = 0.0;
double y = 0.0;
double count = 0.0;
// Sum all members
for (const auto& member : relation->getMembers())
{
if (!map->containsElement(member.getElementId()))
return false;
ConstElementPtr element = map->getElement(member.getElementId());
switch(element->getElementType().getEnum())
{
case ElementType::Node:
{
ConstNodePtr node = std::dynamic_pointer_cast<const Node>(element);
if (!_getCentroidValues(node, map, x, y, count))
return false;
break;
}
case ElementType::Way:
{
ConstWayPtr way = std::dynamic_pointer_cast<const Way>(element);
if (!_getCentroidValues(way, map, x, y, count))
return false;
break;
}
case ElementType::Relation:
{
ConstRelationPtr r = std::dynamic_pointer_cast<const Relation>(element);
if (!_getCentroidValues(r, map, x, y, count))
return false;
break;
}
default:
return false;
}
}
centroid_x += x;
centroid_y += y;
centroid_count += count;
return true;
}

Meters ElementGeometryUtils::_calculateLength(const ConstWayPtr& w, const ConstElementProviderPtr& constProvider)
{
Meters length = 0.0;
const std::vector<long>& node_ids = w->getNodeIds();
for (size_t i = 0; i < node_ids.size() - 1; ++i)
{
// TODO: optimize - we don't really need to convert first, we can just loop through the nodes
// and sum up the distance.
std::shared_ptr<geos::geom::Geometry> geom = ElementToGeometryConverter(constProvider).convertToGeometry(e);
if (geom && geom->isValid())
return geom->getLength();
ConstNodePtr n1 = constProvider->getNode(node_ids[i]);
ConstNodePtr n2 = constProvider->getNode(node_ids[i + 1]);
if (!n1 || !n2)
continue;
length += std::sqrt(std::pow(n2->getX() - n1->getX(), 2.0) + std::pow(n2->getY() - n1->getY(), 2));
}
return 0;
return length;
}

Meters ElementGeometryUtils::_calculateLength(const ConstRelationPtr& r, const ConstElementProviderPtr& constProvider)
{
Meters length = 0.0;
for (const auto& member : r->getMembers())
length += calculateLength(constProvider->getElement(member.getElementId()), constProvider);
return length;
}

}
Loading

0 comments on commit d4276ba

Please sign in to comment.