flopy.utils.gridintersect module

class GridIntersect(mfgrid, rtree=True, local=False)

Bases: object

Class for intersecting shapely geometries with MODFLOW grids.

Notes

• Supports structured and vertex grids. No support for unstructured grids. If grid is layered-unstructured, creating a single layer vertex-grid can be used as a workaround.

• For linestrings and polygons only 2D intersection operations are supported.

• Point intersections can optionally return layer position based on the z-coordinate.

• The STR-tree can be disabled, but this is generally not recommended as some functions will not work without it.

• The STR-tree query is based on the bounding box of the shape or collection, if the bounding box of the shape covers nearly the entire grid, the query won’t be able to limit the search space much, resulting in slower performance. Therefore, it can sometimes be faster to intersect each individual shape in a collection than it is to intersect with the whole collection at once.

filter_query_result(shp, cellids)

Filter a query result to cells that truly intersect a shape.

Used to (further) reduce the query result to cells that intersect with shp. This is primarily used when rtree=False.

Parameters:

• shp (shapely.geometry) – shapely geometry that is prepared and used to filter query result

• cellids (iterable) – iterable of cellids, query result

Returns:

Array of cellids that intersect shp.

Return type:

numpy.ndarray

get_layer_from_z(pts, cellids)

Compute layer indices from point z-values.

Parameters:

• pts (shapely geometry) – Points geometry (single or array-like) with z-values.

• cellids (array_like) – Array of candidate cellids for the points.

Returns:

layer index for each point. Points below/above the grid are returned with index -1.

Return type:

numpy.ndarray

intersect(shp, shapetype=None, sort_by_cellid=True, return_all_intersections=False, contains_centroid=False, min_area_fraction=None, handle_z=False, geo_dataframe=None)

Intersect a shape with a model grid.

Parameters:

• shp (shapely.geometry, geojson object, shapefile.Shape,) – or flopy geometry object

• shapetype (str, optional) – type of shape (i.e. “point”, “linestring”, “polygon” or their multi-variants), used by GeoSpatialUtil if shp is passed as a list of vertices, default is None

• sort_by_cellid (bool) – sort results by cellid, ensures cell with lowest cellid is returned for boundary cases when using vertex methods, default is True

• return_all_intersections (bool, optional) – if True, return multiple intersection results for points or linestrings on grid cell boundaries (e.g. returns 2 intersection results if a point lies on the boundary between two grid cells). The default is False. Only used if shape type is “point” or “linestring”.

• contains_centroid (bool, optional) – if True, only store intersection result if cell centroid is contained within intersection shape, only used if shape type is “polygon”

• min_area_fraction (float, optional) – float defining minimum intersection area threshold, if intersection area is smaller than min_frac_area * cell_area, do not store intersection result, only used if shape type is “polygon”

• handle_z (bool, optional) – Method for handling z-dimension in intersection results for point intersections. Default is False which ignores z-dimension. If True returns the layer index for each point. Points below/above the grid are returned as masked values or pd.NA if returned as dataframe.

• geo_dataframe (bool, optional) – If True, return a geopandas GeoDataFrame, otherwise return a numpy recarray. Default will be set to True in the future; currently, None triggers a deprecation warning and returns a recarray (legacy behavior).

Returns:

Intersection results. Result contains the following columns: - cellid: cellid of intersected grid cell - cellids: legacy column containing (row, col) tuples for structured grids,

or cellids for vertex grids (deprecated, use “cellid” column instead)

• ixshapes or geometry: shapely geometry of the intersection result

And optionally: - layer: layer index of for points with z-dimension when handle_z is True - row: row index of intersected grid cell, for structured grids - col: column index of intersected grid cell, for structured grids - lengths: length of intersection results for linestrings - areas: areas of intersection results for polygons

Return type:

numpy.recarray or geopandas.GeoDataFrame

intersects(shp, shapetype=None, dataframe=None)

Return candidate grid cellids that intersect with shape(s).

Parameters:

• shp (shapely.geometry, geojson geometry, shapefile.shape,) – or flopy geometry object shape to intersect with the grid

• shapetype (str, optional) – type of shape (i.e. “point”, “linestring”, “polygon” or their multi-variants), used by GeoSpatialUtil if shp is passed as a list of vertices, default is None

• dataframe (bool, optional) – If True, return a pandas.DataFrame; otherwise return a numpy recarray. Default will be set to True in the future; currently, None triggers a deprecation warning and returns a recarray (legacy behavior).

Returns:

Grid cell candidates for intersections. Result contains the following columns: - cellid: cellid of candidate grid cell - cellids: legacy column containing (row, col) tuples for structured grids,

or cellids for vertex grids (deprecated, use “cellid” column instead)

And optionally: - row: row index of intersected grid cell, for structured grids - col: column index of intersected grid cell, for structured grids

Return type:

numpy.recarray or pandas.DataFrame

plot_intersection_result(result, plot_grid=True, ax=None, **kwargs)

Plot intersection result.

Parameters:

• result (numpy.rec.recarray or geopandas.GeoDataFrame) – result of intersect()

• plot_grid (bool, optional) – plot model grid, by default True

• ax (matplotlib.Axes, optional) – axes to plot on, by default None which creates a new axis

Returns:

ax – returns axes handle

Return type:

matplotlib.Axes

static plot_linestring(ls, ax=None, cmap=None, **kwargs)

Plot linestrings.

Parameters:

• ls (collection of linestrings) – list, array or GeoSeries containing linestrings

• ax (matplotlib.pyplot.axes, optional) – axes to plot onto, if not provided, creates a new figure

• cmap (str) – matplotlib colormap

• **kwargs – passed to the plot function

Returns:

returns the axes handle

Return type:

matplotlib.pyplot.axes

static plot_point(pts, ax=None, **kwargs)

Plot points.

Parameters:

• pts (collection of points) – array or GeoSeries containing point geometries

• ax (matplotlib.pyplot.axes, optional) – axes to plot onto, if not provided, creates a new figure

• **kwargs – passed to the scatter function

Returns:

returns the axes handle

Return type:

matplotlib.pyplot.axes

static plot_polygon(polys, ax=None, **kwargs)

Plot polygons.

Parameters:

• polys (collection of polygons) – list, array or GeoSeries containing polygons

• ax (matplotlib.pyplot.axes, optional) – axes to plot onto, if not provided, creates a new figure

• **kwargs – passed to the plot function

Returns:

returns the axes handle

Return type:

matplotlib.pyplot.axes

points_to_cellids(pts, handle_z=False, dataframe=True)

Return cellids for points intersecting the grid.

Parameters:

• pts (shapely geometry, geojson geometry, shapefile.Shape, np.ndarray,) – or FloPy geometry object Point(s) to intersect with the grid. array inputs must contain shapely point (or multipoint) geometries.

• handle_z (bool, optional) – Handle z-dimension for points. If True, returns a “layer” column with the computed layer index for each point (NA is returned where the points lie above/below the model grid). Default is False.

• dataframe (bool, optional) – If True, return a pandas.DataFrame; otherwise return a numpy recarray. Default is True.

Returns:

Grid cell indices per point. Result contains the following columns: - cellid: cellid of grid cell containing point - cellids: legacy column containing (row, col) tuples for structured grids,

or cellids for vertex grids (deprecated, use “cellid” column instead)

And optionally: - row: row index of intersected grid cell, for structured grids - col: column index of intersected grid cell, for structured grids - layer: layer index of for points with z-dimension when handle_z is True

Return type:

pandas.DataFrame or numpy.recarray

Notes

Requires rtree=True when initializing GridIntersect.

query_grid(shp, predicate=None)

Perform spatial query on the grid using a shapely geometry.

If no spatial query is possible (rtree=False), returns all grid cells.

Parameters:

• shp (shapely.geometry) – shapely geometry

• predicate (str, optional) – spatial predicate to use for query, default is None. See documentation of self.strtree.query for options.

Returns:

For a single geometry, a 1-D array of cellids. For an array of geometries, a 2xN array where the first row is the input geometry index (shp_id) and the second row contains the matching cellids.

Return type:

numpy.ndarray

parse_shapely_ix_result(collection, ix_result, shptyps=None)

Recursive function for parsing shapely intersection results. Returns a list of shapely shapes matching shptyps.

Parameters:

• collection (list) – state variable for storing result, generally an empty list

• ix_result (shapely.geometry type) – any shapely intersection result

• shptyps (str, list of str, or None, optional) – if None (default), return all types of shapes. if str, return shapes of that type, if list of str, return all types in list

Returns:

list containing shapely geometries of type shptyps

Return type:

list