add support for querying by geometry

[keewis]

No H3 support, and uses shapely.Polygon which might not actually be suitable for this.

[benbovy]

Thanks @keewis!

A couple of comments:

• Do we want to support any shapely.Geometry? This might be quite ambitious regarding (1) support across multiple grid systems (healpy, h3, etc.) and (2) how we actually select grid cells (intersects vs. within the geometry feature).
• If we only support polygons, perhaps easier would be to just pass the coordinates of the polygon vertices and avoid depending on shapely?
• About the name of the methods, what about this pattern?
  • Dataset.dggs.sel_from_points(latitude, longitude)
  • Dataset.dggs.sel_from_geom(geom)
  • Dataset.dggs.sel_from_polygon(interleaved_latlon) or Dataset.dggs.sel_from_polygon(latitude, longitude)

[benbovy]

nit: perhaps useful to leave a small comment about why we need masking.

[benbovy]

About the name of the methods, what about this pattern?

I changed my mind in #20 :). query is consistent with Xvec.

However, I wonder if we shouldn't focus first on methods for extracting cell geometry (#10) or for converting DGGS data cubes into vector data cubes (using xvec). Xvec already implements querying by geometry in a very generic way.

There may be more optimal ways to perform spatial queries on DGGS data than considering it just like any arbitrary vector data? From what I've seen, it looks like DGGS tools are doing things internally that are quite similar to what Xvec does (via shapely / GEOS's rtree).

[keewis]

well, if that's the case then we should probably make use of xvec instead. However, that would mean that we'd have to have query add the boundaries coordinate and index it if it doesn't exist already (or error, but in that case we might as well not provide it)

[benbovy]

However, that would mean that we'd have to have query add the boundaries coordinate and index it if it doesn't exist already (or error, but in that case we might as well not provide it)

Or do not support any spatial indexing in DGGS for now (no Dataset.dggs.query()) and recommend users to do something like below?

# convert to vector data cube
# (automatically set a xvec.GeometryIndex so it is ready for use)
vds = ds.dggs.to_polygons()

# spatial query
vds.xvec.query(...)

That doesn't look like much extra work on the user side.

EDIT: we might as well not provide it as you say.

[keewis]

okay, I agree with not doing this for now until we have evidence that this is possible for every DGGS using just the cell id

[tinaok]

I do want to select the interest area in polygon for my use-case.

[benbovy]

I do want to select the interest area in polygon for my use-case.

Yes, the idea is to support that via Xvec (after converting the DGGS data cube to a polygon vector data cube) and see later if there is a need and/or a more optimal way to support that in xdggs.

[keewis]

it does sound like we're not going to pursue this way of implementing the geometry search any further, so let's close this?

[benbovy]

Yes we can close this I think, and maybe re-open it later if needed.