Tip

For an interactive online version click here:

Convert Wflow staticmaps netcdf to raster files#

In order to inspect or (manually) modify Wflow staticmaps it is convenient to export the maps to a raster format. Here we show how to read the model maps and save to a so-called mapstack (i.e.: a set of raster files with identical grid) using HydroMT.

Load dependencies#

[1]:

import xarray as xr
from os.path import join
import hydromt
from hydromt_wflow import WflowModel

Read wflow staticmaps#

HydroMT provides an easy method to read the model schematization trought the Model API.

[2]:

root = "wflow_piave_subbasin"
mod = WflowModel(root, mode="r")
ds = mod.grid  # here the staticmaps netcdf is loaded
print(ds)

<xarray.Dataset> Size: 1MB
Dimensions:               (latitude: 53, longitude: 58, time: 12, layer: 4)
Coordinates:
  * latitude              (latitude) float64 424B 46.68 46.67 ... 45.83 45.82
  * longitude             (longitude) float64 464B 11.78 11.8 ... 12.72 12.73
  * time                  (time) int32 48B 1 2 3 4 5 6 7 8 9 10 11 12
  * layer                 (layer) int32 16B 0 1 2 3
    spatial_ref           int32 4B 0
Data variables: (12/92)
    x_out                 (latitude, longitude) float64 25kB nan nan ... nan nan
    y_out                 (latitude, longitude) float64 25kB nan nan ... nan nan
    idx_out               (latitude, longitude) int32 12kB -1 -1 -1 ... -1 -1 -1
    wflow_ldd             (latitude, longitude) uint8 3kB 255 255 ... 255 255
    wflow_subcatch        (latitude, longitude) int32 12kB 0 0 0 0 0 ... 0 0 0 0
    wflow_uparea          (latitude, longitude) float32 12kB -9.999e+03 ... -...
    ...                    ...
    TTI                   (latitude, longitude) float32 12kB -999.0 ... -999.0
    TTM                   (latitude, longitude) float32 12kB -999.0 ... -999.0
    WHC                   (latitude, longitude) float32 12kB -999.0 ... -999.0
    G_Cfmax               (latitude, longitude) float32 12kB -999.0 ... -999.0
    G_SIfrac              (latitude, longitude) float32 12kB -999.0 ... -999.0
    G_TT                  (latitude, longitude) float32 12kB -999.0 ... -999.0

Write netcdf to mapstack#

The raster module provides many raster GIS methods throught the raster Dataset accessor. To write a Dataset to a mapstack one line with code is sufficient. We only need to provide the output folder in which all raster files are saved. The default output format is GeoTIFF, but this can be changed with the driver argument. To write to PCRaster map-files it is recommended to have PCRaster python installed.

[3]:

outdir = join(root, "staticmaps")
ds.raster.to_mapstack(outdir)

Now the model files can easily be inspected and modified e.g. QGIS.

NOTE: in QGIS, you can also visualize netcdf files but direct modification is not (yet) possible.

Create staticmaps netcdf files based on mapstack#

If you want to update the staticmaps after modification the maps can be read into a Dataset by hydromt. We recommend the following workflow:

read the original model
read the updated mapstack
change the model root to write the updated model to a new directory
update the staticmaps of the model
write the model

NOTE: We do not read the forcing as it is probably faster to just copy the file instead of loading it into python and writing it back to netcdf.

NOTE: The staticgeoms might be changed because of manual changes in the wflow_river, lakes, reservoir or glacier staticmaps and are therefore not read here. To change these maps we recommend using the hydromt update method to keep the staticgeoms and maps aligned.

[4]:

# read the original model
root = "wflow_piave_subbasin"
mod = WflowModel(root, mode="r")
mod.read_grid()
mod.read_config()

[5]:

# read the updated mapstack
# NOTE: The mapstack does not have to include all staticmaps, only the once that are found will be updated.
# The name of the staticmap should however have to be unchanged.
ds_updated = hydromt.open_mfraster(join(root, "staticmaps", "*.tif"))

[6]:

# change root to a new directory
root_updated = "wflow_piave_subbasin_updated"
mod.set_root(root_updated, mode="w+")

[7]:

# Reinitialize and update the new model staticmaps
mod._grid = xr.Dataset()
mod.set_grid(ds_updated)

[8]:

# write the model to the new directory
mod.write()

Update Wflow staticmaps manually using HydroMT#

The previous steps show you how you can easily save a model staticmaps to a GeoTIFF mapstacks and read it again using hydroMT functions to_mapstack and open_mfraster.

However in order to have a fully ready-to-run wflow model, the mapstacks that we create needs to be processed a little more:

The LAI maps have a third time dimension (other than x and y coordinates)
The c maps have a third soil layer dimension

To update manually a wflow model that is then ready to run, we advise to use the following workflow and functions from the pcraster pcrm module of hydromt_wflow:

read the original model
save the original model staticmaps to formatted PCRaster mapstacks using the write_staticmaps_pcr function (based on to_mapstack with c and LAI pre-processing)
manually update the PCRaster maps (eg using QGIS)
read the updated staticmaps using the read_staticmaps_pcr function (based on open_mfraster with c and LAI post-processing)
Change the model root to write the updated model to a new directory
Write the model

NOTE: We do not read the forcing as it is probably faster to just copy the file instead of loading it into python and writing it back to netcdf.

NOTE: The staticgeoms might be changed because of manual changes in the wflow_river, lakes, reservoir or glacier staticmaps and are therefore not read here. To change these maps we recommend using the hydromt update method to keep the staticgeoms and maps aligned.

Save the staticmaps as a PCRaster mapstack#

[9]:

# read the original model
root = "wflow_piave_subbasin"
mod = WflowModel(root, mode="r+")
mod.read_grid()
mod.read_config()

[10]:

# save the staticmaps to PCRaster mapstacks and update them manually where needed
from hydromt_wflow.pcrm import write_staticmaps_pcr
write_staticmaps_pcr(mod.grid, root=root)

Create a new Wflow model based on the updated PCRaster mapstack#

[11]:

# Import read method for PCRaster files
from hydromt_wflow.pcrm import read_staticmaps_pcr
# create a model instance with the updated staticmaps
root = "wflow_piave_subbasin_updated"
mod = WflowModel(root, mode="w+")
# Read the pcraster staticmaps
staticmaps = read_staticmaps_pcr("wflow_piave_subbasin", crs=4326)
mod.set_grid(staticmaps)
mod.read_config()
# re-generate the basins and rivers staticgeoms, for others (gauges, lakes, reservoirs, glaicers) use the hydromt update method
mod.basins
mod.rivers

[11]:

	geometry	idx	idx_ds	pit	strord
0	LINESTRING (12.73333 46.61667, 12.71667 46.616...	289	346	False	1
1	LINESTRING (12.73333 46.6, 12.71667 46.6)	347	346	False	1
2	LINESTRING (12.7 46.6, 12.71667 46.6)	345	346	False	1
3	LINESTRING (12.28333 46.58333, 12.26667 46.583...	378	492	False	1
4	LINESTRING (12.71667 46.6, 12.73333 46.58333, ...	346	462	False	2
...	...	...	...	...	...
1180	LINESTRING (12.2 45.88333, 12.18333 45.86667)	2809	2866	False	2
1181	LINESTRING (12.16667 45.88333, 12.16667 45.866...	2807	2924	False	1
1182	LINESTRING (12.18333 45.86667, 12.18333 45.85)	2866	2924	False	3
1183	LINESTRING (12.18333 45.85, 12.2 45.85, 12.2 4...	2924	2983	False	5
1184	LINESTRING (12.2 45.83333, 12.2 45.83333)	2983	2983	True	5

1185 rows × 5 columns

[12]:

# write the model to the new directory
mod.write()