fix(mp7particledata): fix get_extent() vertical extent calculation (m…

…odflowpy#2307)

Vertical extent is miscalculated for structured grids when providing a cell index, only the first layer is handled properly. This causes bad results from *ParticleData.to_coords() and .to_prp().

autotest/test_particledata.py

@@ -26,6 +26,7 @@
     ParticleGroupLRCTemplate,
     ParticleGroupNodeTemplate,
 )
+from flopy.modpath.mp7particledata import get_extent
 from flopy.modpath.mp7particlegroup import ParticleGroup
 from flopy.utils.modpathfile import EndpointFile, PathlineFile
 
@@ -47,6 +48,20 @@ def flatten(a):
     ]
 
 
+# test get_extent()
+
+
+def test_get_extent_structured_multilayer():
+    grid = GridCases().structured_small()
+    i, j = 1, 2
+    for k in range(grid.nlay):
+        extent = get_extent(grid, k=k, i=i, j=j)
+        assert extent.minz == grid.botm[k, i, j]
+        assert extent.maxz == (
+            grid.top[i, j] if k == 0 else grid.botm[k - 1, i, j]
+        )
+
+
 # test initializers
 
 

flopy/modpath/mp7particledata.py

@@ -809,7 +809,14 @@ def get_extent(grid, k=None, i=None, j=None, nn=None, localz=False) -> Extent:
     # get cell coords and span in each dimension
     if not (k is None or i is None or j is None):
         verts = grid.get_cell_vertices(i, j)
-        minz, maxz = (0, 1) if localz else (grid.botm[k, i, j], grid.top[i, j])
+        minz, maxz = (
+            (0, 1)
+            if localz
+            else (
+                grid.botm[k, i, j],
+                grid.top[i, j] if k == 0 else grid.botm[k - 1, i, j],
+            )
+        )
     elif nn is not None:
         verts = grid.get_cell_vertices(nn)
         if grid.grid_type == "structured":