move from xugrid noncontiguous to mk.remove_disconnected_regions()

dfm_tools/meshkernel_helpers.py

@@ -161,7 +161,7 @@ def geographic_to_meshkernel_projection(is_geographic:bool) -> meshkernel.Projec
     return projection
 
 
-def meshkernel_to_UgridDataset(mk:meshkernel.MeshKernel, crs:(int,str) = None, remove_noncontiguous:bool = False) -> xu.UgridDataset:
+def meshkernel_to_UgridDataset(mk:meshkernel.MeshKernel, crs:(int,str) = None) -> xu.UgridDataset:
     """
     Convert a meshkernel object to a UgridDataset, including a variable with the crs (used by dflowfm to distinguish spherical/cartesian networks).
     The UgridDataset enables bathymetry interpolation and writing to netfile.
@@ -172,8 +172,6 @@ def meshkernel_to_UgridDataset(mk:meshkernel.MeshKernel, crs:(int,str) = None, r
         DESCRIPTION.
     crs : (int,str), optional
         DESCRIPTION. The default is None.
-    remove_noncontiguous : bool, optional
-        DESCRIPTION. The default is False.
 
     Returns
     -------
@@ -188,23 +186,6 @@ def meshkernel_to_UgridDataset(mk:meshkernel.MeshKernel, crs:(int,str) = None, r
 
     xu_grid = xu.Ugrid2d.from_meshkernel(mesh2d_grid)
 
-    #remove non-contiguous grid parts
-    def xugrid_remove_noncontiguous(grid):
-        #based on https://deltares.github.io/xugrid/examples/connectivity.html#connected-components
-        #uses https://docs.scipy.org/doc/scipy/reference/sparse.csgraph.html
-        #TODO: maybe replace with meshkernel?
-        uda = xu.UgridDataArray(
-            xr.DataArray(np.ones(grid.node_face_connectivity.shape[0]), dims=["face"]), grid
-        )
-        labels = uda.ugrid.connected_components()
-        counts = labels.groupby(labels).count()
-        most_frequent_label = counts["group"][np.argmax(counts.data)].item() #find largest contiguous part
-        labels = labels.where(labels == most_frequent_label, drop=True)
-        grid = labels.grid
-        return grid
-    if remove_noncontiguous:
-        xu_grid = xugrid_remove_noncontiguous(xu_grid)
-
     #convert 0-based to 1-based indices for connectivity variables like face_node_connectivity
     xu_grid_ds = xu_grid.to_dataset()
     xu_grid_ds = xr.decode_cf(xu_grid_ds) #decode_cf is essential since it replaces fillvalues with nans

tests/examples_workinprogress/workinprogress_meshkernel_creategrid.py

@@ -15,7 +15,7 @@
 import datetime as dt
 
 #general settings
-lon_min,lon_max = -6,2
+lon_min,lon_max = -5,2
 lat_min,lat_max = 48.5,51.2
 lon_res,lat_res = 0.2,0.2
 projection = meshkernel.ProjectionType.SPHERICAL #True for spherical grids
@@ -106,11 +106,25 @@
 dfmt.plot_coastlines(ax=ax, res='h', min_area=100)
 ctx.add_basemap(ax=ax, crs=crs, attribution=False)
 
+
+"""
+remove disconnected parts of grid (from non-contiguous to contiguous)
+"""
+
+mk.mesh2d_remove_disconnected_regions()
+
+mesh2d_contiguous = mk.mesh2d_get()
+fig, ax = plt.subplots(figsize=figsize)
+mesh2d_contiguous.plot_edges(ax,linewidth=0.8)
+dfmt.plot_coastlines(ax=ax, res='h', min_area=100)
+ctx.add_basemap(ax=ax, crs=crs, attribution=False)
+
+
 """
 convert meshkernel grid to xugrid, interp bathymetry, plot and save to *_net.nc
 """
 
-xu_grid_uds = dfmt.meshkernel_to_UgridDataset(mk, remove_noncontiguous=True, crs=crs) #TODO: put remove_noncontiguous in meshkernel?: https://github.com/Deltares/MeshKernelPy/issues/44
+xu_grid_uds = dfmt.meshkernel_to_UgridDataset(mk, crs=crs)
 
 fig, ax = plt.subplots(figsize=figsize)
 xu_grid_uds.grid.plot(ax=ax,linewidth=0.8) #TODO: maybe make uds instead of ds (but then bathy interpolation goes wrong)