Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Add centroid, boundary and convex_hull function (returning new Geography) #20

Merged
merged 14 commits into from
Aug 19, 2024
Merged

Add centroid, boundary and convex_hull function (returning new Geography) #20

Show file tree
Hide file tree
Changes from 4 commits
Commits
Show all changes
14 commits
Select commit Hold shift + click to select a range
3d7dd92
basic support for convex_hull
jorisvandenbossche Mar 8, 2023
496d6aa
Merge remote-tracking branch 'upstream/main' into concave_hull
jorisvandenbossche Mar 8, 2023
4f65fdf
add tests
jorisvandenbossche Mar 8, 2023
0f5fb03
cleanup
jorisvandenbossche Mar 9, 2023
1cd9376
from_geog method + return Polygon
jorisvandenbossche Mar 10, 2023
82a7e71
Update src/accessors-geog.cpp
jorisvandenbossche Mar 18, 2023
3331305
Merge remote-tracking branch 'upstream/main' into concave_hull
jorisvandenbossche Mar 18, 2023
5a5d714
update tests
jorisvandenbossche Mar 18, 2023
8f5a025
linting
jorisvandenbossche Mar 18, 2023
11b0fbb
add centroid and boundary
jorisvandenbossche Mar 18, 2023
91ce98f
add type annotations
jorisvandenbossche Mar 18, 2023
b899e07
Merge remote-tracking branch 'upstream/main' into concave_hull
jorisvandenbossche Aug 19, 2024
c368271
add to api docs
jorisvandenbossche Aug 19, 2024
22e1439
Merge remote-tracking branch 'upstream/main' into concave_hull
jorisvandenbossche Aug 19, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
1 change: 1 addition & 0 deletions CMakeLists.txt
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -52,6 +52,7 @@ endif()

add_library(spherely MODULE
src/geography.cpp
src/accessors-geog.cpp
src/predicates.cpp
src/spherely.cpp
)
Expand Down
29 changes: 29 additions & 0 deletions src/accessors-geog.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1,29 @@
#include <s2geography.h>

#include "geography.hpp"
#include "pybind11.hpp"

namespace py = pybind11;
namespace s2geog = s2geography;
using namespace spherely;

PyObjectGeography convex_hull(PyObjectGeography a) {
const auto& a_ptr = a.as_geog_ptr()->geog();
auto s2_obj = s2geog::s2_convex_hull(a_ptr);
auto geog_ptr = std::make_unique<Geography>(std::move(s2_obj));
Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

So currently, this function returns instances of the base class Geography. If I change this to std::make_unique<Polygon>, then it actually returns a Polygon python subclass.

For the convex_hull case, I think I can be sure that this will always return a Polygon. But in general, do we need a method on s2geography::Geography to know which subtype it actually is (a PolygonGeography) in this case), so we can assert / known to which subclass to downcast on our side?

Copy link
Owner

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think the subtype could be inferred from both num_shapes() and dimension() methods, except maybe for empty geographies. That said, a specific method might be indeed convenient and slightly more efficient as dimension() iterates over each shape in "multi" subclasses.

I'm not sure if we really need it, though. Do you have some cases in mind where this would need to be inferred at runtime?

For the convex hull case, s2geog::s2_convex_hull returns a s2geog::Geography pointer but it should probably return a s2geog::PolygonGeography pointer since this is what S2ConvexHullAggregator::Finalize() returns anyway.

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I don't know if this is the same for s2, but with GEOS there are certainly geometry-returning functions that can return varying geometry types depending on the input (eg an intersection can basically return any geometry type)

Copy link
Owner

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I don't know if this is the same for s2, but with GEOS there are certainly geometry-returning functions that can return varying geometry types depending on the input (eg an intersection can basically return any geometry type)

If you know the type of the input(s) can you always determine the type of the output(s), or is there some geometry-returning functions where the output geometry type depends on some internal processing logic?

In the 1st case, we could use spherely::Geography::geog_type() returned from the inputs. Although I agree this might be worth of having such method available in s2geography, which would help getting rid of spherely::Geography wrapper classes that makes things quite confused.

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

If you know the type of the input(s) can you always determine the type of the output(s), or is there some geometry-returning functions where the output geometry type depends on some internal processing logic?

Yes, AFAIK that is the case. For example consider intersection of two polygons: if they overlap, you can get a Polygon or MultiPolygon, if they touch along a segment you get a LineString (or MultiLineString), and if they touch you get a Point (or MultiPoint), and if they both overlap in some area and touch in another area, you get a GeometryCollection.

auto res_object = PyObjectGeography::as_py_object(std::move(geog_ptr));
return static_cast<PyObjectGeography&>(res_object);
}

void init_accessors(py::module& m) {
m.def("convex_hull", py::vectorize(&convex_hull), py::arg("a"),
R"pbdoc(
Computes the convex hull of each geography.

Parameters
----------
a : :py:class:`Geography` or array_like
Geography object

)pbdoc");
}
4 changes: 4 additions & 0 deletions src/pybind11.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -90,6 +90,10 @@ class PyObjectGeography : public py::object {
return py::cast(std::move(geog_ptr));
}

// static PyObjectGeography as_py_object(std::unique_ptr<T> geog_ptr) {
// return static_cast<PyObjectGeography>(py::cast(std::move(geog_ptr)));
// }
Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I need to clean this up, but I thought to have as_py_object actually return PyObjectGeography instead of py::object, since it's currently not yet used, and in practice we need PyObjectGeography in the vectorized functions. But didn't get this to work ..

Copy link
Owner

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Something like this should work:

template <class T, std::enable_if<...>>
static PyObjectGeography as_py_object(std::unique_ptr<T> geog_ptr) {
    auto pyobj = py::cast(std::move(geog_ptr));
    auto pyobj_geog = static_cast<PyObjectGeography&>(pyobj);
    return std::move(pyobj_geog);
}

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thanks, will give that a try. Are you OK with changing the existing as_py_object to do that? (I am not sure if we have a use case of the current version that returns py::object)

Copy link
Owner

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yes sure, it is also more consistent to have a static method returning an instance of the same type IMO. (from_geog could be a better name than as_py_object).


// Just check whether the object is a Geography
//
bool is_geog_ptr() const { return check_type(false); }
Expand Down
2 changes: 2 additions & 0 deletions src/spherely.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -7,6 +7,7 @@ namespace py = pybind11;

void init_geography(py::module&);
void init_predicates(py::module&);
void init_accessors(py::module&);

PYBIND11_MODULE(spherely, m) {
m.doc() = R"pbdoc(
Expand All @@ -19,6 +20,7 @@ PYBIND11_MODULE(spherely, m) {

init_geography(m);
init_predicates(m);
init_accessors(m);

#ifdef VERSION_INFO
m.attr("__version__") = MACRO_STRINGIFY(VERSION_INFO);
Expand Down
32 changes: 32 additions & 0 deletions tests/test_accessors.py
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1,32 @@
import numpy as np

import spherely

import pytest


@pytest.mark.parametrize(
"geog, expected",
[
# (spherely.Point(0, 0), spherely.Point(0, 0)),
(
spherely.LineString([(0, 0), (0, 2), (2, 2)]),
spherely.Polygon([(0, 0), (0, 2), (2, 2)])
),
# (
# spherely.Polygon([(0, 0), (0, 2), (2, 2), (1, 0.5)]),
# spherely.Polygon([(0, 0), (0, 2), (2, 2)])
# ),
],
)
def test_convex_hull(geog, expected) -> None:
# test array + scalar
actual = spherely.convex_hull(geog)
assert isinstance(actual, spherely.Geography)
assert spherely.equals(actual, expected)

actual = spherely.convex_hull([geog])
assert isinstance(actual, np.ndarray)
actual = actual[0]
assert isinstance(actual, spherely.Geography)
assert spherely.equals(actual, expected)