Proposal to consolidate geospatial inputs

[mnfienen]

Proposal for consolidating geospatial object input handling

There are a handful of methods and functions that take various geospatial geometries as input (e.g. points, lines, or polygons). Currently, these are handled in various ways, typically not taking advantage of the shapely.geometry types of Point, Line, or Polygon. Most functions take a list of points that describe a line or a polygon and a separate argument defining which type it represents.

The goal of this issue is to consolidate all this handling with a helper function that can intercept geometries and convert them to whatever type the underlying function/method needs.

The approach would be to allow users to supply geometry as one of the following:

• GeoJSON
• shapefile
• lists of vertices
• shapely.geomtry objects
• shapely.geometry collections
• custom flopy classes for geometry from geometry.py
• lists of shapely.geometry objects (e.g. the geometry column from a geodataframe)

A key issue will be handling the dependency on shapely. If we can somehow isolate that dependency in the helper that intercepts the geometry we won't need to import it elsewhere. Unless, of course, we just allow shapely to be a base dependency of flopy throughout.

The functions and methods listed below take geospatial information as inputs. These would all need an interface to the helper function to parse inputs and make necessary translations.

Parent File	Function/Method	Current options
triangle.py	add_polygon	list of (x,y) points
gridgen.py	add_refinement_features	list of shapes or shapefile name - need to verify types
gridgen.py	add_active_domain	list of shapes or shapefile name - need to verify
gridgen.py	intersect	list of shapes or shapefile name - need to verify
geometry.py	get_polygon_area, get_polygon_centroid, is_clockwise	list of vertices
geometry.py	Point, LineString, Polygon types	pyshp types - we should expand to allow shapely.geometry types. One way to optionally support shapely types would be to inherit from them if shapely can be imported. That would allow us to keep current attributes and methods while adding shapely functionality. Alternatively, the non-shapely attributes and methods on these could be implemented elsewhere.
gridintersect.py	_intersect_linestring_shapely, _intersect_linestring_structured, _intersect_point_shapely, _intersect_point_structured, _intersect_polygon_shapely, _intersect_polygon_structured,	These all seem to use shapely.geometry objects - perhaps can consolidate into one function call that accepts any of these types
rasters.py	sample_polygon, _intersection,_sample_rio_dataset, crop,	shapely.geometry or GeoJSON
rasters.py	_point_in_polygon,	List of vertices
rasters.py	sample_point,	x and y coordinates
shapefile_utils.py	recarray2shp,	geometry types defined in geometry.py. recarray2shp is used in mfsfr2.py, utils.py and modpathfile.py
shapefile_utils.py	shp2recarray,	shp2recarray reads a shapefile and casts the feature geometries to the types in geometry.py
grid module	nothing at this time	in my opinion it would be useful to have a property method on modelgrid instances that would fetch (and cache) cell geometries

[jdhughes-dev]

@jlarsen-usgs @dbrakenhoff and thoughts on this suggestion?

[jlarsen-usgs]

@jdhughes-usgs @dbrakenhoff I think this would be a nice addition to FloPy and would make working from spatial objects much easier. I can start putting together a set of utilities for this proposal.

[dbrakenhoff]

@mnfienen @jlarsen-usgs This sounds like a really good idea! There's quite a bit of code in GridIntersect related to extracting coordinates from all the different flavors of shapely geometries. If we implement your suggestion, I'm sure that could simplified. Additionally, I like the idea of users being able to pass all manner of geometries to flopy methods. I'll gladly take a look at anything you propose!

[jlarsen-usgs]

@dbrakenhoff and @mnfienen
I started putting together a base to work from that allows for easy capture (within flopy) and conversion from one geospatial input type to another. This class is fairly simple and takes advantage of the __geo_interface__ (geojson standard) method that is attached to objects in pyshp and in shapely. The class is not yet complete (I haven't added custom flopy geometry classes or methods to capture and convert geometry collections/feature collections yet), but I thought I'd share it to begin the feedback process.

My thought on shapefile, geometry collection, and feature collection conversion is that there will be a thin wrapper class (much like util3d wraps around util2d) to handle bulk conversions.

My idea for handling custom flopy geometries is to add the geojson standard __geo_interface__ property and a _from_geojson() method for easy conversion.

I've attached a copy of the module here if you'd like to take a look. It's been tested with points, polygons, polygons with holes, multipolygons, and multipolygons with holes so far.

geospatial_utils.zip