diff --git a/src/geography.cpp b/src/geography.cpp index 6432d9b..92d58cf 100644 --- a/src/geography.cpp +++ b/src/geography.cpp @@ -1,13 +1,5 @@ #include "geography.hpp" -#include -#include -#include -#include -#include -#include -#include - #include #include #include @@ -16,6 +8,14 @@ #include #include +#include +#include +#include +#include +#include +#include +#include + #include "pybind11.hpp" #include "pybind11/stl.h" @@ -38,13 +38,11 @@ py::detail::type_info *PyObjectGeography::geography_tinfo = nullptr; // It is not fully supported with Pybind11 (Point is non-default constructible) // https://github.com/pybind/pybind11/issues/4108 // -S2Point to_s2point(const std::pair& vertex) { +S2Point to_s2point(const std::pair &vertex) { return S2LatLng::FromDegrees(vertex.first, vertex.second).ToPoint(); } -S2Point to_s2point(const Point* vertex) { - return vertex->s2point(); -} +S2Point to_s2point(const Point *vertex) { return vertex->s2point(); } /* ** Helper to create Geography object wrappers. @@ -54,12 +52,11 @@ S2Point to_s2point(const Point* vertex) { ** @tparam S The S2Geometry type */ template -std::unique_ptr make_geography(S&& s2_obj) { +std::unique_ptr make_geography(S &&s2_obj) { S2GeographyPtr s2geog_ptr = std::make_unique(std::forward(s2_obj)); return std::make_unique(std::move(s2geog_ptr)); } - class PointFactory { public: static std::unique_ptr FromLatLonDegrees(double lat_degrees, @@ -69,29 +66,30 @@ class PointFactory { return make_geography(S2Point(latlng)); } - //TODO: from LatLonRadians + // TODO: from LatLonRadians }; template -static std::unique_ptr create_linestring(const std::vector& coords) { +static std::unique_ptr create_linestring( + const std::vector &coords) { std::vector pts(coords.size()); - std::transform( - coords.begin(), coords.end(), pts.begin(), - [](const V& vertex) { return to_s2point(vertex); }); + std::transform(coords.begin(), coords.end(), pts.begin(), + [](const V &vertex) { return to_s2point(vertex); }); auto polyline_ptr = std::make_unique(pts); - return make_geography(std::move(polyline_ptr)); + return make_geography( + std::move(polyline_ptr)); } template -static std::unique_ptr create_polygon(const std::vector& shell) { +static std::unique_ptr create_polygon( + const std::vector &shell) { std::vector shell_pts(shell.size()); - std::transform( - shell.begin(), shell.end(), shell_pts.begin(), - [](const V& vertex) { return to_s2point(vertex); }); + std::transform(shell.begin(), shell.end(), shell_pts.begin(), + [](const V &vertex) { return to_s2point(vertex); }); auto shell_loop_ptr = std::make_unique(); // TODO: maybe add an option to skip validity checks @@ -118,7 +116,8 @@ static std::unique_ptr create_polygon(const std::vector& s polygon_ptr->set_s2debug_override(S2Debug::DISABLE); polygon_ptr->InitOriented(std::move(loops)); - return make_geography(std::move(polygon_ptr)); + return make_geography( + std::move(polygon_ptr)); } /* @@ -132,7 +131,7 @@ py::array_t num_shapes(const py::array_t geographies) { py::buffer_info result_buf = result.request(); int *rptr = static_cast(result_buf.ptr); - for (size_t i = 0; i < buf.size; i++) { + for (py::ssize_t i = 0; i < buf.size; i++) { auto geog_ptr = (*geographies.data(i)).as_geog_ptr(); rptr[i] = geog_ptr->num_shapes(); } @@ -157,9 +156,7 @@ py::array_t create(py::array_t xs, double *yptr = static_cast(ybuf.ptr); py::object *rptr = static_cast(rbuf.ptr); - size_t size = static_cast(xbuf.shape[0]); - - for (size_t i = 0; i < xbuf.shape[0]; i++) { + for (py::ssize_t i = 0; i < xbuf.shape[0]; i++) { auto point_ptr = PointFactory::FromLatLonDegrees(xptr[i], yptr[i]); // rptr[i] = PyObjectGeography::as_py_object(std::move(point_ptr)); rptr[i] = py::cast(std::move(point_ptr)); @@ -283,10 +280,9 @@ void init_geography(py::module &m) { pylinestring.def(py::init(&create_linestring>), py::arg("coordinates")); - pylinestring.def(py::init(&create_linestring), + pylinestring.def(py::init(&create_linestring), py::arg("coordinates")); - auto pypolygon = py::class_(m, "Polygon", R"pbdoc( A geography type representing an area that is enclosed by a linear ring. @@ -301,8 +297,9 @@ void init_geography(py::module &m) { )pbdoc"); - pypolygon.def(py::init(&create_polygon>), py::arg("shell")); - pypolygon.def(py::init(&create_polygon), py::arg("shell")); + pypolygon.def(py::init(&create_polygon>), + py::arg("shell")); + pypolygon.def(py::init(&create_polygon), py::arg("shell")); // Temp test diff --git a/src/geography.hpp b/src/geography.hpp index 05627b1..f4bb60f 100644 --- a/src/geography.hpp +++ b/src/geography.hpp @@ -14,7 +14,12 @@ using S2GeographyIndexPtr = std::unique_ptr; /* ** The registered Geography types */ -enum class GeographyType : std::int8_t { None = -1, Point, LineString, Polygon }; +enum class GeographyType : std::int8_t { + None = -1, + Point, + LineString, + Polygon +}; /* ** Thin wrapper around s2geography::Geography. @@ -81,7 +86,8 @@ class Point : public Geography { } inline const S2Point& s2point() const { - const auto& points = static_cast(geog()).Points(); + const auto& points = + static_cast(geog()).Points(); // TODO: does not work for empty point geography return points[0]; } diff --git a/src/predicates.cpp b/src/predicates.cpp index 5e7c5d0..4f8275b 100644 --- a/src/predicates.cpp +++ b/src/predicates.cpp @@ -1,6 +1,8 @@ #include #include +#include + #include "geography.hpp" #include "pybind11.hpp" @@ -8,48 +10,32 @@ namespace py = pybind11; namespace s2geog = s2geography; using namespace spherely; -bool intersects(PyObjectGeography a, PyObjectGeography b) { - const auto& a_index = a.as_geog_ptr()->geog_index(); - const auto& b_index = b.as_geog_ptr()->geog_index(); - - S2BooleanOperation::Options options; - return s2geog::s2_intersects(a_index, b_index, options); -} - -bool equals(PyObjectGeography a, PyObjectGeography b) { - const auto& a_index = a.as_geog_ptr()->geog_index(); - const auto& b_index = b.as_geog_ptr()->geog_index(); - - S2BooleanOperation::Options options; - return s2geog::s2_equals(a_index, b_index, options); -} - -bool contains(PyObjectGeography a, PyObjectGeography b) { - const auto& a_index = a.as_geog_ptr()->geog_index(); - const auto& b_index = b.as_geog_ptr()->geog_index(); - - S2BooleanOperation::Options options; - return s2geog::s2_contains(a_index, b_index, options); -} - -bool within(PyObjectGeography a, PyObjectGeography b) { - const auto& a_index = a.as_geog_ptr()->geog_index(); - const auto& b_index = b.as_geog_ptr()->geog_index(); - - S2BooleanOperation::Options options; - return s2geog::s2_contains(b_index, a_index, options); -} +/* +** Functor for predicate bindings. +*/ +class Predicate { +public: + using FuncType = std::function; + + template + Predicate(F&& func) : m_func(std::forward(func)) {} + + bool operator()(PyObjectGeography a, PyObjectGeography b) const { + const auto& a_index = a.as_geog_ptr()->geog_index(); + const auto& b_index = b.as_geog_ptr()->geog_index(); + return m_func(a_index, b_index, m_options); + } -bool disjoint(PyObjectGeography a, PyObjectGeography b) { - const auto& a_index = a.as_geog_ptr()->geog_index(); - const auto& b_index = b.as_geog_ptr()->geog_index(); - - S2BooleanOperation::Options options; - return !s2geog::s2_intersects(a_index, b_index, options); -} +private: + FuncType m_func; + S2BooleanOperation::Options m_options; +}; void init_predicates(py::module& m) { - m.def("intersects", py::vectorize(&intersects), py::arg("a"), py::arg("b"), + m.def("intersects", py::vectorize(Predicate(s2geog::s2_intersects)), + py::arg("a"), py::arg("b"), R"pbdoc( Returns True if A and B share any portion of space. @@ -62,7 +48,8 @@ void init_predicates(py::module& m) { )pbdoc"); - m.def("equals", py::vectorize(&equals), py::arg("a"), py::arg("b"), + m.def("equals", py::vectorize(Predicate(s2geog::s2_equals)), py::arg("a"), + py::arg("b"), R"pbdoc( Returns True if A and B are spatially equal. @@ -76,7 +63,8 @@ void init_predicates(py::module& m) { )pbdoc"); - m.def("contains", py::vectorize(&contains), py::arg("a"), py::arg("b"), + m.def("contains", py::vectorize(Predicate(s2geog::s2_contains)), + py::arg("a"), py::arg("b"), R"pbdoc( Returns True if B is completely inside A. @@ -86,9 +74,16 @@ void init_predicates(py::module& m) { Geography object(s) )pbdoc"); - - m.def("within", py::vectorize(&within), py::arg("a"), py::arg("b"), - R"pbdoc( + + m.def( + "within", + py::vectorize(Predicate([](const s2geog::ShapeIndexGeography& a_index, + const s2geog::ShapeIndexGeography& b_index, + const S2BooleanOperation::Options& options) { + return s2geog::s2_contains(b_index, a_index, options); + })), + py::arg("a"), py::arg("b"), + R"pbdoc( Returns True if A is completely inside B. Parameters @@ -98,8 +93,15 @@ void init_predicates(py::module& m) { )pbdoc"); - m.def("disjoint", py::vectorize(&disjoint), py::arg("a"), py::arg("b"), - R"pbdoc( + m.def( + "disjoint", + py::vectorize(Predicate([](const s2geog::ShapeIndexGeography& a_index, + const s2geog::ShapeIndexGeography& b_index, + const S2BooleanOperation::Options& options) { + return !s2geog::s2_intersects(a_index, b_index, options); + })), + py::arg("a"), py::arg("b"), + R"pbdoc( Returns True if A boundaries and interior does not intersect at all with those of B. @@ -109,5 +111,4 @@ void init_predicates(py::module& m) { Geography object(s) )pbdoc"); - - } +} diff --git a/src/pybind11.hpp b/src/pybind11.hpp index 0b39d81..f48b2f3 100644 --- a/src/pybind11.hpp +++ b/src/pybind11.hpp @@ -142,7 +142,8 @@ struct npy_format_descriptor { }; // Override signature type hint for vectorized Geography arguments -template struct handle_type_name> { +template +struct handle_type_name> { static constexpr auto name = _("Geography | array_like"); };