cleaned up meshkernel_delete_withshp and testcase (#566)

Deltares · Oct 5, 2023 · 53205a1 · 53205a1
1 parent ac4691c
commit 53205a1
Show file tree

Hide file tree

Showing 8 changed files with 38 additions and 13 deletions.
diff --git a/README.md b/README.md
@@ -21,4 +21,3 @@ A Python package for pre- and postprocessing D-FlowFM model input and output fil
 - install with ``pip install dfm_tools -U`` (or [installation guide](https://deltares.github.io/dfm_tools/installation))
 - [online documentation](https://deltares.github.io/dfm_tools) with installation guide, contributing guide, tutorials/examples, API reference and a convenient search box.
 - Bug or feature request? Create a [GitHub issue](https://github.com/Deltares/dfm_tools/issues)
-
diff --git a/dfm_tools/meshkernel_helpers.py b/dfm_tools/meshkernel_helpers.py
@@ -48,17 +48,20 @@ def meshkernel_delete_withcoastlines(mk:meshkernel.meshkernel.MeshKernel, res:st
     coastlines_gdf = get_coastlines_gdb(bbox=bbox, res=res, min_area=min_area, crs=crs)
 
     meshkernel_delete_withgdf(mk, coastlines_gdf)
-
-def meshkernel_delete_withshp(mk:meshkernel.meshkernel.MeshKernel, coastlines_shp, min_area: float=0, crs:(int,str) = None):
+
+
+def meshkernel_delete_withshp(mk:meshkernel.meshkernel.MeshKernel, coastlines_shp:str, crs:(int,str) = None):
     """
     Delete parts of mesh that are inside the shapefile polygon.
 
     Parameters
     ----------
     mk : meshkernel.meshkernel.MeshKernel
         DESCRIPTION.
-    coastlines_shp : path to the shp file
-        DESCRIPTION.
+    coastlines_shp : str
+        Path to the shp file.
+    crs : (int,str), optional
+        DESCRIPTION. The default is None.
 
     Returns
     -------
@@ -67,16 +70,15 @@ def meshkernel_delete_withshp(mk:meshkernel.meshkernel.MeshKernel, coastlines_sh
     """
 
     mesh_bnds = mk.mesh2d_get_mesh_boundaries_as_polygons()
-
     bbox = (mesh_bnds.x_coordinates.min(), mesh_bnds.y_coordinates.min(), mesh_bnds.x_coordinates.max(), mesh_bnds.y_coordinates.max())
-
-    coastlines_gdb = gpd.read_file(coastlines_shp, include_fields= ['polygons'], bbox=bbox)
+    coastlines_gdb = gpd.read_file(coastlines_shp, bbox=bbox)
 
     if crs:
         coastlines_gdb = coastlines_gdb.to_crs(crs)
 
     meshkernel_delete_withgdf(mk, coastlines_gdb)
-
+
+
 def meshkernel_delete_withgdf(mk:meshkernel.meshkernel.MeshKernel, coastlines_gdf:gpd.GeoDataFrame):
     """
     Delete parts of mesh that are inside the polygons/Linestrings in a GeoDataFrame.

diff --git a/tests/shp/test.cpg b/tests/shp/test.cpg
diff --git a/tests/shp/test.dbf b/tests/shp/test.dbf
diff --git a/tests/shp/test.prj b/tests/shp/test.prj
diff --git a/tests/shp/test.shp b/tests/shp/test.shp
diff --git a/tests/shp/test.shx b/tests/shp/test.shx
diff --git a/tests/test_meshkernel_helpers.py b/tests/test_meshkernel_helpers.py
@@ -12,6 +12,9 @@
 import meshkernel
 import xarray as xr
 import numpy as np
+import geopandas as gpd
+from shapely.geometry import Polygon
+import glob
 
 
 @pytest.mark.unittest
@@ -71,8 +74,26 @@ def test_meshkernel_delete_withcoastlines():
 
     assert len(mk.mesh2d_get().face_nodes) == 17364
 
+
 @pytest.mark.unittest
 def test_meshkernel_delete_withshp():
+    # write shapefile from coords
+    file_shp = 'mk_delete_test.shp'
+    points_x = np.array([-68.40156631636155, -68.36143523088661, -68.28392176442131, -68.26413109213229,  
+                         -68.20915700244058, -68.1965129618115, -68.20860726154366, -68.199811407193,    
+                         -68.23059689742034, -68.23389534280184, -68.26303161033846, -68.29436684146273, 
+                         -68.27512591007063, -68.3064611411949, -68.37462901241263, -68.41146165250606, 
+    					 -68.41311087519682, -68.40156631636155])
+    points_y = np.array([12.301796772360385, 12.303445995051137, 12.243524237287176, 12.233628901142668, 
+                         12.225382787688911, 12.206141856296814, 12.184152220420131, 12.115984349202414, 
+                         12.080800931799722, 12.028025805695682, 12.033523214664854, 12.113785385614745, 
+                         12.144021134945184, 12.20119418822456,  12.220984860513575, 12.226482269482746, 
+                         12.286404027246707, 12.301796772360385])
+
+    geom = Polygon(zip(points_x, points_y))
+    gdf = gpd.GeoDataFrame(geometry=[geom], crs='EPSG:4326')
+    gdf.to_file(file_shp)
+
     #generate basegrid
     lon_min, lon_max, lat_min, lat_max = -68.45, -68.1, 12, 12.35
     dxy = 0.005
@@ -81,9 +102,14 @@ def test_meshkernel_delete_withshp():
     assert len(mk.mesh2d_get().face_nodes) == 20732
 
     # remove cells with a shapefile
-    dfmt.meshkernel_delete_withshp(mk=mk, coastlines_shp='./shp/test.shp')
+    dfmt.meshkernel_delete_withshp(mk=mk, coastlines_shp=file_shp)
 
-    assert len(mk.mesh2d_get().face_nodes) == 17180
+    assert len(mk.mesh2d_get().face_nodes) == 17272
+
+    # delete shapefile
+    shp_list = glob.glob(file_shp.replace('.shp','.*'))
+    [os.remove(x) for x in shp_list]
+
 
 @pytest.mark.unittest
 def test_meshkernel_delete_withgdf():