Convert dfs2 to raster

[adrienmichez]

Hi,

Would you have some code to share about the conversion from dfs2 to raster grids?

I tried to do so using rioxarray from this example but with no success

Isn't there a easier way to convert mikeio grid to raster file? I would like to find a python solution similar as the GUI tool in mikezero toolbox (GIS/Mike2Grd).

Best
Adrien

[Snowthe]

I support this - might be a good new feature if a to_raster() or to_geotiff() method could be added to the mikeio.Dataset class

Until then, I wrote a little function for myself that converts a mikeio.Dataset or DataArray (one timestep and layer) to a geotiff file - maybe that can serve as inspiration:

def dfs_to_geotiff(ds, out_fp, item=None, time=None, layer=None):
    """
    Writes a mikeio Dataset/DataArray to geotiff. Requires specification of 
    item, time and layer in case a whole dataset with multiple dimensions is 
    passed.
    The output geotiff is simply a 2D raster (without multiple bands/layers).
    Uses rioxarray.

    Parameters
    ----------
    ds : mikeio.Dataset or mikeio.DataArray
        to be written to geotiff
    out_fp : str
        File path of the geotiff output
    item : int or str, optional
        If default=None: Use first item
        If int: Use specified item number
        If string: Use specified item name
    time : datetime or int, optional
        If default=None: Use first timestep
        If int: Use specified timestep index
        If datetime: Use timestep closest to specified time
    layer : int, optional
        If default=None: Use uppermost layer
        if int: Use specified layer
        positive values: 0 gives the lowest layer, <no_layers-1> gives the top layer
        negative values: -1 gives the top layer, -<no_layers> gives the lowest layer
    Returns
    -------
    None.
    """
    # only import rioxarray here if actually used
    import rioxarray
    if item==None: #get right item
        item = 0
    elif (isinstance(item, int)) | (isinstance(item, str)):
        pass
    else:
        sys.exit('ERROR: item has to be int, str or None')
    if time==None: #get time index
        ti = 0
    elif isinstance(time, int):
        ti = time
    elif isinstance(time, datetime):
        ti = _get_tsi_from_mikeio(ds, time)[0]
    else:
        sys.exit('ERROR: time has to be int, datetime or None')
    if layer==None: #get right layer
        li = -1 #uppermost layer
    elif isinstance(layer, int):
        pass
    else:
        sys.exit('ERROR: layer has to be int or None')
    # obtain data to be written to geotiff
    if ('z' in ds.dims) & ('time' in ds.dims): #if 3D data with both time and z dim
        out_data = ds[item].to_xarray()[ti, li, :, :]
    elif ('z' in ds.dims) & len(ds.dims==3): #if 3D data with z dim, but without time (can happen if timestep is "pre-selected")
        out_data = ds[item].to_xarray()[li, :, :]
    elif len(ds.dims==3): #assuming that else, dimensions are x, y, time
        out_data = ds[item].to_xarray()[ti, :, :]
    elif len(ds.dims==2): #if 2D data, but without time (can happen if timestep is "pre-selected")
        out_data = ds[item].to_xarray()[:, :]
    # write to geotiff
    out_data.rio.to_raster(out_fp)