From a0e9219407b7be208d38f4f902cd2b5b96da1351 Mon Sep 17 00:00:00 2001
From: wpbonelli <wbonelli@ucar.edu>
Date: Thu, 12 Sep 2024 11:15:03 -0400
Subject: [PATCH] fix(mp7particledata): fix get_extent() vertical extent
 calculation (#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    | 15 +++++++++++++++
 flopy/modpath/mp7particledata.py |  9 ++++++++-
 2 files changed, 23 insertions(+), 1 deletion(-)

diff --git a/autotest/test_particledata.py b/autotest/test_particledata.py
index f3032cb6e2..f5363f957d 100644
--- a/autotest/test_particledata.py
+++ b/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
 
 
diff --git a/flopy/modpath/mp7particledata.py b/flopy/modpath/mp7particledata.py
index f905219398..472ff19151 100644
--- a/flopy/modpath/mp7particledata.py
+++ b/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":