Tip

For an interactive online version click here: Binder badge

Add (local) emission data#

Once you have a DELWAQ model, you may want to update your model in order to add new emission data, add sample locations, use different hydrological forcing data, create and run different scenarios etc.

With HydroMT, you can easily read your model and update one or several components of your model using the update function of the command line interface (CLI). Here are the steps and some examples on how to update and add local emission data to your model.

All lines in this notebook which starts with ! are executed from the command line. Within the notebook environment the logging messages are shown after completion. You can also copy these lines and paste them in your shell to get more feedback.

Import packages#

In this notebook, we will use some functions of HydroMT to prepare configuration and catalog files and plot the new emission data of the updated model. Here are the libraries to import to realize these steps.

[1]:
import numpy as np
import geopandas as gpd
[2]:
# for plotting
import matplotlib.pyplot as plt
import cartopy.io.img_tiles as cimgt
import cartopy.crs as ccrs
proj = ccrs.PlateCarree() # plot projection
[3]:
# import hydromt
import hydromt
from hydromt_delwaq import DelwaqModel, DemissionModel
[4]:
# setup logging
from  hydromt.log import setuplog
logger = setuplog("update_model_emission_local", log_level=10)
2024-06-21 02:53:46,926 - update_model_emission_local - log - INFO - HydroMT version: 0.10.0

HydroMT CLI update interface#

Using the HydroMT update API, we can update one or several components of an already existing DELWAQ model. Let’s get an overview of the available options:

[5]:
# Print the options available from the update command
! hydromt update --help
Usage: hydromt update [OPTIONS] MODEL MODEL_ROOT

  Update a specific component of a model.

  Set an output directory to copy the edited model to a new folder, otherwise
  maps are overwritten.

  Example usage: --------------

  Update (overwrite!) landuse-landcover based maps in a Wflow model: hydromt
  update wflow /path/to/model_root -c setup_lulcmaps --opt lulc_fn=vito -d
  /path/to/data_catalog.yml -v

  Update Wflow model components outlined in an .yml configuration file and
  write the model to a directory: hydromt update wflow /path/to/model_root  -o
  /path/to/model_out  -i /path/to/wflow_config.yml  -d
  /path/to/data_catalog.yml -v

Options:
  -o, --model-out DIRECTORY  Output model folder. Maps in MODEL_ROOT are
                             overwritten if left empty.
  -i, --config PATH          Path to hydroMT configuration file, for the model
                             specific implementation.
  -c, --components TEXT      Model methods from configuration file to run
  --opt TEXT                 Method specific keyword arguments, see the method
                             documentation of the specific model for more
                             information about the arguments.
  -d, --data TEXT            Path to local yaml data catalog file OR name of
                             predefined data catalog.
  --dd, --deltares-data      Flag: Shortcut to add the "deltares_data" catalog
  --fo, --force-overwrite    Flag: If provided overwrite existing model files
  --cache                    Flag: If provided cache tiled rasterdatasets
  -q, --quiet                Decrease verbosity.
  -v, --verbose              Increase verbosity.
  --help                     Show this message and exit.

Adding local data to the model#

Using HydroMT, it is rather easy to add additional emission data using local data. Here we will see an example where we add a emission factor map from a local (dummy!) source stored in examples_data/emission_factor.gpkg.

Compared to global data, that can be added directly to a delwaq model, we first need to add our local data the the HydroMT Data Catalog using a local yaml data library. The full steps and infomation to add data the HydroMT data library is available in the docs.

Here we will see the steps for our (dummy) emission factor source. And first see what the data looks like (you can also download and open it in QGIS instead of loading it with python).

[6]:
#Let's read and quickly plot the data
emission_data = 'examples_data/emission_factor.gpkg'
#Open with geopandas and plot
gdf = gpd.read_file(emission_data)
#Content
print("Content of the local emission_factor.gpkg file: ")
gdf
Content of the local emission_factor.gpkg file:
[6]:
OBJECTID ISO_3DIGIT NAME EF geometry
0 19 AUT Austria 1.500 MULTIPOLYGON (((15.01626 49.01902, 15.03378 49...
1 22 BEL Belgium 2.260 MULTIPOLYGON (((4.94299 51.44159, 4.94078 51.4...
2 47 CHE Switzerland 1.710 MULTIPOLYGON (((8.57693 47.81010, 8.59278 47.8...
3 65 DEU Germany 1.390 MULTIPOLYGON (((8.70214 47.71345, 8.71030 47.7...
4 82 FRA France 2.175 MULTIPOLYGON (((9.25835 41.34536, 9.25813 41.3...
5 135 LIE Liechtenstein 1.605 MULTIPOLYGON (((9.56471 47.17015, 9.57504 47.1...
6 139 LUX Luxembourg 2.260 MULTIPOLYGON (((6.02926 50.17731, 6.03060 50.1...
7 171 NLD Netherlands 1.410 MULTIPOLYGON (((3.51562 51.40709, 3.54988 51.4...
8 117 ITA Italy 2.175 MULTIPOLYGON (((12.53929 35.52953, 12.54022 35...
[7]:
#coordinate system
print(f"Coordinate system of the data: {gdf.crs}")
#Quick plot
gdf.plot()
Coordinate system of the data: EPSG:4326
[7]:
<Axes: >
../_images/_examples_adding_local_emission_14_2.png

We can see that our data is a vector file with some European countries in EPSG 4326, and that it contains our emission factors in the EF column. We now have enough information to add it to the HydroMT data sources by preparing a local_emission_sources.yml.

[8]:
# Dictionary with all the required information on our emission_factor.gpkg file
data_dict = {
    'EF_local': {                      # user defined internal name of the local data source
        'path': 'examples_data/emission_factor.gpkg', # path to the local data
        'data_type': 'GeoDataFrame',   # hydroMT DataCatalog type 'GeoDataFrame' for vector file
        'driver': 'vector',            # driver to read the file 'vector' for vector file
        'crs': 4326,                   # optional here but mentioned as en example
        'rename': {
            'NAME': 'COUNTRY',         # dummy here but can be used to rename some of the data columns
        },
        'unit_mult': {
            'EF': 1.0,                 # dummy here (EF*1.0) but can be used to convert units (should be kg/d/EV)
        },
        'unit_add':{
            'EF': 0.0,                 # dummy here (EF+0.0) but can be used to convert units (should be kg/d/EV)
        },
        'meta': {                      # additional information on the file (license, download link, DOI, author...)
            'category': 'socio econpmic',
            'source_url': 'https://github.com/Deltares/hydromt_delwaq/tree/main/examples/examples_data/emission_factor.gpkg',
            'notes': 'Dummy emission factor data',
        },
    },
}

# Convert the dict to hydroMT yaml library format and save the file
fn_yml = 'local_emission_sources.yml'
data_catalog = hydromt.DataCatalog(logger=logger)
data_catalog.from_dict(data_dict)
data_catalog.to_yml(fn_yml)

#Read the saved yml
with open(fn_yml, 'r') as f:
    txt = f.read()
print(txt)
meta:
  root: /home/runner/work/hydromt_delwaq/hydromt_delwaq/docs/_examples
EF_local:
  data_type: GeoDataFrame
  path: examples_data/emission_factor.gpkg
  driver: vector
  rename:
    NAME: COUNTRY
  unit_mult:
    EF: 1.0
  unit_add:
    EF: 0.0
  meta:
    category: socio econpmic
    source_url: https://github.com/Deltares/hydromt_delwaq/tree/main/examples/examples_data/emission_factor.gpkg
    notes: Dummy emission factor data
  crs: 4326

Our local data is now ready to be processed by HydroMT !

Preparing the configuration file#

As our new EF_local data is a vector data, we will add it to our model using the [setup_emission_vector] component of HydroMT Delwaq.

Let’s prepare a HydroMT configuration file with our options for preparing the EF_local grid. All available options are available in the docs(setup_emission_vector).

[9]:
# Dictionnary with all the components and options we want to update
emission_vector_options = {
    'setup_emission_vector': {
        'emission_fn': 'EF_local',
        'col2raster': 'EF',
        'rasterize_method': 'value',
    },
}

# Save it to a hydroMT ini file
fn_ini = "delwaq_update_emission_local.yml"
hydromt.config.configwrite(fn_ini, emission_vector_options)

# Open the file and visualize the content
with open(fn_ini, 'r') as f:
    txt = f.read()
print(txt)
setup_emission_vector:
  emission_fn: EF_local
  col2raster: EF
  rasterize_method: value

Some explanations about the option we chose here:

  • emission_fn: name of the emission factor data source in HydroMT Data Catalog. The one we choose when creating the local_emission_sources.yml.

  • col2raster: name of the column in the vector file that contains the emission factors values to rasterize to the model grid.

  • rasterize_method: method to rasterize the vector either value to rasterize the values in the ‘col2raster’ column, or fraction to rasterize the fraction of the model grid cell that is covered by the vector shapes (eg fraction of agricultural areas).

Updating the model with the local data#

[10]:
# NOTE: copy this line (without !) to your shell for more direct feedback
! hydromt update demission EM_piave -o ./EM_piave_extended_local -i delwaq_update_emission_local.yml -d local_emission_sources.yml --fo -vv
2024-06-21 02:53:52,433 - update - log - DEBUG - Writing log messages to new file /home/runner/work/hydromt_delwaq/hydromt_delwaq/docs/_examples/EM_piave_extended_local/hydromt.log.
2024-06-21 02:53:52,433 - update - log - INFO - HydroMT version: 0.10.0
2024-06-21 02:53:52,433 - update - main - INFO - Updating demission model at /home/runner/work/hydromt_delwaq/hydromt_delwaq/docs/_examples/EM_piave (r).
2024-06-21 02:53:52,433 - update - main - INFO - Output dir: /home/runner/work/hydromt_delwaq/hydromt_delwaq/docs/_examples/EM_piave_extended_local
2024-06-21 02:53:52,433 - update - main - INFO - User settings:
2024-06-21 02:53:52,450 - update - data_catalog - INFO - Parsing data catalog from local_emission_sources.yml
2024-06-21 02:53:52,451 - update - log - DEBUG - Appending log messages to file /home/runner/work/hydromt_delwaq/hydromt_delwaq/docs/_examples/EM_piave/hydromt.log.
2024-06-21 02:53:52,451 - update - model_api - INFO - Initializing demission model from hydromt_delwaq (v0.2.2.dev0).
2024-06-21 02:53:52,451 - update - model_api - ERROR - Model config file not found at /home/runner/work/hydromt_delwaq/hydromt_delwaq/docs/_examples/EM_piave/demission.inp
2024-06-21 02:53:52,452 - update - model_api - DEBUG - Reading model file monareas.
2024-06-21 02:53:52,480 - update - model_api - DEBUG - Reading model file staticdata.
2024-06-21 02:53:52,525 - update - model_api - DEBUG - Reading model file basins.
nodata value missing for /home/runner/work/hydromt_delwaq/hydromt_delwaq/docs/_examples/EM_piave/hydromodel/rivarea.tif
2024-06-21 02:53:52,653 - update - model_api - DEBUG - Reading model file staticdata.
2024-06-21 02:53:52,661 - update - model_grid - WARNING - Replacing grid map: river
2024-06-21 02:53:52,662 - update - model_grid - WARNING - Replacing grid map: slope
2024-06-21 02:53:52,662 - update - model_grid - WARNING - Replacing grid map: streamorder
2024-06-21 02:53:52,663 - update - model_grid - WARNING - Replacing grid map: SoilThickness
2024-06-21 02:53:52,664 - update - model_grid - WARNING - Replacing grid map: porosity
2024-06-21 02:53:52,664 - update - model_grid - WARNING - Replacing grid map: ptiddown
2024-06-21 02:53:52,665 - update - model_grid - WARNING - Replacing grid map: monareas
2024-06-21 02:53:52,666 - update - model_grid - WARNING - Replacing grid map: ghs_pop_2015
2024-06-21 02:53:52,667 - update - model_grid - WARNING - Replacing grid map: gdp_world
2024-06-21 02:53:52,667 - update - demission - INFO - Model read
2024-06-21 02:53:52,668 - update - log - DEBUG - Appending log messages to file /home/runner/work/hydromt_delwaq/hydromt_delwaq/docs/_examples/EM_piave_extended_local/hydromt.log.
2024-06-21 02:53:52,675 - update - model_api - INFO - setup_emission_vector.emission_fn: EF_local
2024-06-21 02:53:52,675 - update - model_api - INFO - setup_emission_vector.col2raster: EF
2024-06-21 02:53:52,675 - update - model_api - INFO - setup_emission_vector.rasterize_method: value
2024-06-21 02:53:52,675 - update - delwaq - INFO - Preparing 'EF_local' map.
2024-06-21 02:53:52,677 - update - geodataframe - INFO - Reading EF_local vector data from /home/runner/work/hydromt_delwaq/hydromt_delwaq/docs/_examples/examples_data/emission_factor.gpkg
2024-06-21 02:53:52,896 - update - geodataframe - DEBUG - Clip intersects [11.775, 45.808, 12.742, 46.692] (EPSG:4326)
2024-06-21 02:53:52,898 - update - geodataframe - DEBUG - Convert units for 2 columns.
2024-06-21 02:53:52,980 - update - demission - INFO - Write model data to /home/runner/work/hydromt_delwaq/hydromt_delwaq/docs/_examples/EM_piave_extended_local
2024-06-21 02:53:53,007 - update - delwaq - INFO - Writing staticmap files.
2024-06-21 02:53:53,145 - update - model_api - DEBUG - Writing file geoms/monareas.geojson
2024-06-21 02:53:53,156 - update - model_api - DEBUG - Writing file geoms/basins.geojson
2024-06-21 02:53:53,162 - update - delwaq - INFO - Writing model config to file.
2024-06-21 02:53:53,163 - update - delwaq - INFO - Writing hydromap files.
nodata value missing for /home/runner/work/hydromt_delwaq/hydromt_delwaq/docs/_examples/EM_piave_extended_local/hydromodel/rivarea.tif
2024-06-21 02:53:53,280 - update - demission - DEBUG - No forcing data found, skip writing.

The example above means the following: run hydromt with:

  • update demission: i.e. update a demission model

  • EM_piave: original model folder

  • -o ./EM_piave_extended_local: output updated model folder

  • -i delwaq_update_emission_local.yml: hydroMT configuration file containing the components and options to update

  • -d local_emission_sources.yml: hydroMT local data library file containing local data sources and their attributes.

  • -v: give some extra verbosity (2 * v) to display feedback on screen. Now debug messages are provided.

Visualization of the new emission map#

We can now plot our newly created emission factor map.

[11]:
# Load the original and updated model with hydromt
mod = DemissionModel(root='EM_piave_extended_local', mode='r')
[12]:
# Edit the lines below to change the emission map and its colormap
emissionmap = 'EF_local'
colormap = 'Reds'

#Load the emission map
da = mod.grid[emissionmap].raster.mask_nodata()
da.attrs.update(long_name=emissionmap, units='-')

#Plot
# we assume the model maps are in the geographic CRS EPSG:4326
proj = ccrs.PlateCarree()
# adjust zoomlevel and figure size to your basis size & aspect
zoom_level = 10
figsize=(10, 8)

# initialize image with geoaxes
fig = plt.figure(figsize=figsize)
ax = fig.add_subplot(projection=proj)
extent = np.array(da.raster.box.buffer(0.02).total_bounds)[[0, 2, 1, 3]]
ax.set_extent(extent, crs=proj)

# add sat background image
ax.add_image(cimgt.QuadtreeTiles(), zoom_level, alpha=0.5)

# add the country shape
country_shape = 'examples_data/emission_factor.gpkg'
#Open with geopandas and plot
gdf = gpd.read_file(country_shape)
gdf.boundary.plot(transform=proj, ax=ax, label="Country")

## plot emission map
cmap = plt.cm.get_cmap(colormap)
kwargs = dict(cmap=cmap)
# plot 'normal' elevation
da.plot(transform=proj, ax=ax, zorder=1, cbar_kwargs=dict(aspect=30, shrink=.8), **kwargs)

ax.xaxis.set_visible(True)
ax.yaxis.set_visible(True)
ax.set_ylabel(f"latitude [degree north]")
ax.set_xlabel(f"longitude [degree east]")
_ = ax.set_title(f"Delwaq emission map")
legend = ax.legend(
    handles=[*ax.get_legend_handles_labels()[0]],
    title="Legend",
    loc='lower right',
    frameon=True,
    framealpha=0.7,
    edgecolor='k',
    facecolor='white'
)
../_images/_examples_adding_local_emission_28_0.png

Our Delwaq model is completely located in Italy so our EF_local map here only has the same value for every grid cell.

If you are using Binder, feel free to edit this notebook to try out with your own local data (you need to first upload your data in Binder using the upload button).

Note that you can also download the models you created in Binder by first zipping them using the lines below and downloading the zipped model.

[13]:
# Lines to zip a model folder
model_folder = 'EM_piave_extended_local'

#Zipping
import shutil
shutil.make_archive(model_folder, 'zip', model_folder)
[13]:
'/home/runner/work/hydromt_delwaq/hydromt_delwaq/docs/_examples/EM_piave_extended_local.zip'