Skip to content

bug: GridIntersect's intersect method may throw IndexError when LineString starts outside model grid #2049

Description

@vincentpost

This error can (not will) occur when a LineString starts outside the model area for model grids that have a rotated grid or an x,y offset. I ran into it when using FloPy version 3.3.5.

It is hard to reproduce but as far as I can see the bug is caused after the LineString gets clipped (in my case in _intersect_linestring_structured, line 821 of grid_intersect.py: lineclip = shp.intersection(pl)). When _get_nodes_intersecting_linestring is called next, the first point of the clipped line is transformed to real-world coordinates using transform(...) and then intersect is called on line 1009: (i, j) = self.intersect(shapely_geo.Point(x0[0], y0[0])).cellids[0] to determine the cellid of the point. The IndexError results when the intersect method returns an empty recarray. This can happen because when the coordinates of the first point of the line are back-transformed to local coordinates (within _intersect_point_structured line 725 in gridintersect.py: rx, ry = transform(... ), the local x coordinate rx is not 0 but (in my case) -1.4551915228366852e-11. As a result, the point falls outside the bounds of the grid in local coordinates and the returned recarray is empty, causing the IndexError.

A possible solution might be to use math.isclose() within find_position_in_array in gridintersect.py. But before changing the code I thought I'd throw it out here for discussion.

Metadata

Metadata

Assignees

No one assigned

    Labels

    Type

    No type

    Projects

    No projects

    Milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions