From 2647ca0939e07d5fac76696c09fdb522a428d88b Mon Sep 17 00:00:00 2001 From: Daniel Baston Date: Wed, 19 Jun 2024 13:17:07 -0400 Subject: [PATCH] Curved geometries: Implement getLength() --- include/geos/algorithm/CircularArcs.h | 1 - include/geos/geom/CircularArc.h | 45 +++++++++++++-- include/geos/geom/CircularString.h | 6 +- src/geom/CircularString.cpp | 18 ++++++ tests/unit/capi/GEOSLengthTest.cpp | 4 +- tests/unit/geom/CircularArcTest.cpp | 76 ++++++++++++++++++++++++++ tests/unit/geom/CircularStringTest.cpp | 4 +- tests/unit/geom/CompoundCurveTest.cpp | 4 +- tests/unit/geom/CurvePolygonTest.cpp | 2 +- tests/unit/geom/MultiCurveTest.cpp | 4 +- tests/unit/geom/MultiSurfaceTest.cpp | 3 +- 11 files changed, 148 insertions(+), 19 deletions(-) create mode 100644 tests/unit/geom/CircularArcTest.cpp diff --git a/include/geos/algorithm/CircularArcs.h b/include/geos/algorithm/CircularArcs.h index 7e3b9b7bb..54f0a9b7d 100644 --- a/include/geos/algorithm/CircularArcs.h +++ b/include/geos/algorithm/CircularArcs.h @@ -23,7 +23,6 @@ namespace algorithm { class GEOS_DLL CircularArcs { public: - static constexpr double PI = 3.14159265358979323846; /// Return the circle center of an arc defined by three points static geom::CoordinateXY getCenter(const geom::CoordinateXY& p0, const geom::CoordinateXY& p1, diff --git a/include/geos/geom/CircularArc.h b/include/geos/geom/CircularArc.h index 758733738..69ca66746 100644 --- a/include/geos/geom/CircularArc.h +++ b/include/geos/geom/CircularArc.h @@ -47,7 +47,7 @@ class GEOS_DLL CircularArc { int orientation() const { if (!m_orientation_known) { - m_orientation = algorithm::Orientation::index(center(), p0, p1); + m_orientation = algorithm::Orientation::index(p0, p1, p2); m_orientation_known = true; } return m_orientation; @@ -71,6 +71,43 @@ class GEOS_DLL CircularArc { return m_radius; } + bool isLinear() const { + return std::isnan(radius()); + } + + bool isCircle() const { + return p0.equals(p2); + } + + double angle() const { + if (isCircle()) { + return 2*MATH_PI; + } + + auto t0 = theta0(); + auto t2 = theta2(); + + if (orientation() == algorithm::Orientation::COUNTERCLOCKWISE) { + std::swap(t0, t2); + } + + if (t0 < t2) { + t0 += 2*MATH_PI; + } + + auto diff = t0-t2; + + return diff; + } + + double length() const { + if (isLinear()) { + return p0.distance(p2); + } + + return angle()*radius(); + } + double theta0() const { return std::atan2(p0.y - center().y, p0.x - center().x); } @@ -123,11 +160,11 @@ class GEOS_DLL CircularArc { t2 -= t0; theta -= t0; - if (t2 < 0){ - t2 += 2*algorithm::CircularArcs::PI; + if (t2 < 0) { + t2 += 2*MATH_PI; } if (theta < 0) { - theta += 2*algorithm::CircularArcs::PI; + theta += 2*MATH_PI; } return theta >= t2; diff --git a/include/geos/geom/CircularString.h b/include/geos/geom/CircularString.h index 283cb636c..fd43bc2b9 100644 --- a/include/geos/geom/CircularString.h +++ b/include/geos/geom/CircularString.h @@ -15,7 +15,6 @@ #pragma once #include -#include namespace geos { namespace geom { @@ -35,10 +34,7 @@ class GEOS_DLL CircularString : public SimpleCurve { GeometryTypeId getGeometryTypeId() const override; - double getLength() const override - { - throw util::UnsupportedOperationException("Cannot calculate length of CircularString"); - } + double getLength() const override; bool hasCurvedComponents() const override { diff --git a/src/geom/CircularString.cpp b/src/geom/CircularString.cpp index f5bb375c8..7a758f7a1 100644 --- a/src/geom/CircularString.cpp +++ b/src/geom/CircularString.cpp @@ -12,6 +12,7 @@ * **********************************************************************/ +#include #include #include #include @@ -49,6 +50,23 @@ CircularString::getGeometryTypeId() const return GEOS_CIRCULARSTRING; } +double +CircularString::getLength() const +{ + if (isEmpty()) { + return 0; + } + + const CoordinateSequence& coords = *getCoordinatesRO(); + + double tot = 0; + for (std::size_t i = 2; i < coords.size(); i += 2) { + auto len = CircularArc(coords[i-2], coords[i-1], coords[i]).length(); + tot += len; + } + return tot; +} + CircularString* CircularString::reverseImpl() const { diff --git a/tests/unit/capi/GEOSLengthTest.cpp b/tests/unit/capi/GEOSLengthTest.cpp index 66c2872d1..7be71de71 100644 --- a/tests/unit/capi/GEOSLengthTest.cpp +++ b/tests/unit/capi/GEOSLengthTest.cpp @@ -1,6 +1,7 @@ #include // geos #include +#include #include "capi_test_utils.h" @@ -67,7 +68,8 @@ void object::test<4>() double length = -1; int ret = GEOSLength(input_, &length); - ensure_equals("error raised on curved geometry", ret, 0); + ensure_equals(ret, 1); + ensure_equals(length, geos::MATH_PI); } } // namespace tut diff --git a/tests/unit/geom/CircularArcTest.cpp b/tests/unit/geom/CircularArcTest.cpp new file mode 100644 index 000000000..00f7bbcbd --- /dev/null +++ b/tests/unit/geom/CircularArcTest.cpp @@ -0,0 +1,76 @@ +#include + +#include +#include +#include + +using geos::geom::CoordinateXY; +using geos::geom::CircularArc; +using geos::MATH_PI; + +namespace tut { + +struct test_circulararc_data { + + const double eps = 1e-8; + + void checkAngle(const CoordinateXY& p0, const CoordinateXY& p1, const CoordinateXY& p2, double expected) { + CircularArc arc(p0, p1, p2); + ensure_equals(p0.toString() + " / " + p1.toString() + " / " + p2.toString(), arc.angle(), expected, eps); + + CircularArc rev(p2, p1, p0); + ensure_equals(p2.toString() + " / " + p1.toString() + " / " + p0.toString(), rev.angle(), expected, eps); + } + + void checkLength(const CoordinateXY& p0, const CoordinateXY& p1, const CoordinateXY& p2, double expected) { + CircularArc arc(p0, p1, p2); + ensure_equals(p0.toString() + " / " + p1.toString() + " / " + p2.toString(), arc.length(), expected, eps); + + CircularArc rev(p2, p1, p0); + ensure_equals(p2.toString() + " / " + p1.toString() + " / " + p0.toString(), rev.length(), expected, eps); + } +}; + +using group = test_group; +using object = group::object; + +group test_circulararc_group("geos::geom::CircularArc"); + +// test angle() on unit circle +template<> +template<> +void object::test<1>() +{ + auto x = std::sqrt(2.0)/2; + + // full circle + checkAngle({-1, 0}, {1, 0}, {-1, 0}, 2*MATH_PI); + + // check half-circles + checkAngle({-1, 0}, {0, 1}, {1, 0}, MATH_PI); // top + checkAngle({-1, 0}, {0, -1}, {1, 0}, MATH_PI); // bottom + checkAngle({0, -1}, {-1, 0}, {0, 1}, MATH_PI); // left + checkAngle({0, -1}, {1, 0}, {0, 1}, MATH_PI); // right + + // check quadrants + checkAngle({-1, 0}, {-x, x}, {0, 1}, MATH_PI/2); // upper left + checkAngle({0, 1}, {x, x}, {1, 0}, MATH_PI/2); // upper right + checkAngle({0, -1}, {x, -x}, {1, 0}, MATH_PI/2); // lower right + checkAngle({0, -1}, {-x, -x}, {-1, 0}, MATH_PI/2); // lower left + + // 3/4 + checkAngle({-x, x}, {0, -1}, {x, x}, 1.5*MATH_PI); // mouth up + checkAngle({-x, -x}, {0, 1}, {x, -x}, 1.5*MATH_PI); // mouth down + checkAngle({-x, x}, {1, 0}, {-x, -x}, 1.5*MATH_PI); // mouth left + checkAngle({x, -x}, {-1, 0}, {x, x}, 1.5*MATH_PI); // mouth right +} + +// test length() +template<> +template<> +void object::test<2>() +{ + checkLength({1.6, 0.4}, {1.6, 0.5}, {1.7, 1}, 0.6122445326877711); +} + +} diff --git a/tests/unit/geom/CircularStringTest.cpp b/tests/unit/geom/CircularStringTest.cpp index 8e1eaf4f3..64bbe6524 100644 --- a/tests/unit/geom/CircularStringTest.cpp +++ b/tests/unit/geom/CircularStringTest.cpp @@ -58,7 +58,7 @@ void object::test<1>() ensure(cs->getCoordinate() == nullptr); ensure_equals(cs->getArea(), 0); - ensure_THROW(cs_->getLength(), geos::util::UnsupportedOperationException); + ensure_equals(cs->getLength(), 0); } // Basic Geometry API @@ -74,7 +74,7 @@ void object::test<2>() // Geometry size functions ensure("isEmpty", !cs_->isEmpty()); ensure_equals("getArea", cs_->getArea(), 0); - ensure_THROW(cs_->getLength(), geos::util::UnsupportedOperationException); + ensure_equals("getLength", cs_->getLength(), 2*geos::MATH_PI); ensure_equals("getNumGeometries", cs_->getNumGeometries(), 1u); ensure_equals("getNumPoints", cs_->getNumPoints(), 5u); geos::geom::Envelope expected(0, 4, -1, 1); diff --git a/tests/unit/geom/CompoundCurveTest.cpp b/tests/unit/geom/CompoundCurveTest.cpp index ae232f5ca..a430842fa 100644 --- a/tests/unit/geom/CompoundCurveTest.cpp +++ b/tests/unit/geom/CompoundCurveTest.cpp @@ -65,7 +65,7 @@ void object::test<1>() ensure("getCoordinate", cc->getCoordinate() == nullptr); ensure_equals("getArea", cc->getArea(), 0); - ensure_THROW(cc_->getLength(), geos::util::UnsupportedOperationException); + ensure_equals("getLength", cc->getLength(), 0); } // Basic Geometry API @@ -81,7 +81,7 @@ void object::test<2>() // Geometry size functions ensure("isEmpty", !cc_->isEmpty()); ensure_equals("getArea", cc_->getArea(), 0); - ensure_THROW(cc_->getLength(), geos::util::UnsupportedOperationException); + ensure_equals("getLength", cc_->getLength(), geos::MATH_PI + 2); ensure_equals("getNumGeometries", cc_->getNumGeometries(), 1u); ensure_equals("getNumCurves", cc_->getNumCurves(), 2u); ensure_equals("getNumPoints", cc_->getNumPoints(), 5u); // FIXME should this be 5 or 4? diff --git a/tests/unit/geom/CurvePolygonTest.cpp b/tests/unit/geom/CurvePolygonTest.cpp index 0025235d9..9cafc6516 100644 --- a/tests/unit/geom/CurvePolygonTest.cpp +++ b/tests/unit/geom/CurvePolygonTest.cpp @@ -91,7 +91,7 @@ void object::test<2>() // Geometry size functions ensure("isEmpty", !cp_->isEmpty()); ensure_THROW(cp_->getArea(), geos::util::UnsupportedOperationException); - ensure_THROW(cp_->getLength(), geos::util::UnsupportedOperationException); + ensure_equals("getLength", cp_->getLength(), 19.236489581872586, 1e-8); ensure_equals("getNumGeometries", cp_->getNumGeometries(), 1u); ensure_equals("getNumPoints", cp_->getNumPoints(), 14u); ensure_equals("getNumInteriorRing", cp_->getNumInteriorRing(), 1u); diff --git a/tests/unit/geom/MultiCurveTest.cpp b/tests/unit/geom/MultiCurveTest.cpp index 1fec0ebc2..ec7844c8b 100644 --- a/tests/unit/geom/MultiCurveTest.cpp +++ b/tests/unit/geom/MultiCurveTest.cpp @@ -81,7 +81,7 @@ void object::test<1>() ensure("getCoordinate", mc->getCoordinate() == nullptr); ensure_equals("getArea", mc->getArea(), 0); - ensure_THROW(mc_->getLength(), geos::util::UnsupportedOperationException); + ensure_equals("getLength", mc->getLength(), 0); } // Basic Geometry API @@ -97,7 +97,7 @@ void object::test<2>() // Geometry size functions ensure("isEmpty", !mc_->isEmpty()); ensure_equals("getArea", mc_->getArea(), 0); - ensure_THROW(mc_->getLength(), geos::util::UnsupportedOperationException); + ensure_equals("getLength", mc_->getLength(), 22.064916706618778, 1e-8); ensure_equals("getNumGeometries", mc_->getNumGeometries(), 3u); ensure_equals("getNumPoints", mc_->getNumPoints(), 16u); ensure(!mc_->getEnvelopeInternal()->isNull()); diff --git a/tests/unit/geom/MultiSurfaceTest.cpp b/tests/unit/geom/MultiSurfaceTest.cpp index 08c51f9e5..751010c23 100644 --- a/tests/unit/geom/MultiSurfaceTest.cpp +++ b/tests/unit/geom/MultiSurfaceTest.cpp @@ -7,6 +7,7 @@ #include #include #include +#include using geos::geom::CoordinateXY; @@ -85,7 +86,7 @@ void object::test<2>() // Geometry size functions ensure("isEmpty", !ms_->isEmpty()); ensure_THROW(ms_->getArea(), geos::util::UnsupportedOperationException); - ensure_THROW(ms_->getLength(), geos::util::UnsupportedOperationException); + ensure_equals("getLength", ms_->getLength(), 10.283185307179586); ensure_equals("getNumGeometries", ms_->getNumGeometries(), 2u); ensure_equals("getNumPoints", ms_->getNumPoints(), 10u); ensure(!ms_->getEnvelopeInternal()->isNull());