Content:
Installation
Preparing model input parameters and dataset
Runing DRYP model
Post processing model outputs
Preprocessing datasets
This section will help you to create dataset required as model parameters for runing DRYP. This section uses simple raster operations to store and save variables (scripts can be accessed at https://github.com/AndresQuichimbo/DRYPv2.0.1.git.). More complex operations can be performed using any GIS software (e.g ArcMap, QGis) o python libraries such as rasterio.
For this section we will use ready-to-use datasets created for running the regional model. Raw data can be found in different repositories such as Hydrosheds, SoilGrids, etc.
The following steps will be carried out in order to produce the set of model paramter datasets for a subbasin of the HAD region.
Delineate a new subbasin from the regional model domain
Clip extent of model parameters to reduce simulation times (reduce model domain)
Create a flow direction map compatible with DRYP-landlab formart
Create a river network
Create a constant head boundary condition from water table depth
The preprocessing (also postprocessing) module is under developement therefore the list of functions available are still limited. All available functions can be found in the following link: https://github.com/AndresQuichimbo/DRYPv2.0.1.git.
[ ]:
# Import libraries from local repository
import sys
sys.path.append('C:/Users/Edisson/Documents/GitHub/DRYPv2.0.1') # CHANGE THIS TO YOUR MODEL PATH
import cuwalid.dryp.components.DRYP_watershed as ppbasin
import cuwalid.tools.DRYP_pptools as pptools
import cuwalid.tools.DRYP_rrtools as rrtools
NOTE: Before we start, please download and unzip all files required for this training by cliking in the following link: https://drive.google.com/drive/folders/1Lf9fdOjjZO87sU1PCsUbKTNjt5IXr9DD?usp=drive_link
The file is organized as bellow:
[ ]:
print("Local_directory\n|__training\n |__basin\n | |__datasets\n | | |__csv\n | | |__shp\n | |__model\n | | |__inputs\n | |__outputs\n |__regional")
WARNING: Once the training files are downloaded change the paths below to access them during this training material.
[ ]:
training_general_path = "D:/HAD/training/"
regional_path = "D:/HAD/training/regional/"
basin_path = "D:/HAD/training/basin/"
1. Delineate a watershed (getting a basin mask)
For basin delineation, you will need to provide a Digital Elevation Map (DEM), the flow direction, and the basin outlet as csv file
First, specify path of parameter files and create maps for each one of them.
[ ]:
fname_surface = regional_path + "/model/inputs/HAD_DEM_utm_mm.asc"
fname_flowdir = regional_path + "/model/inputs/HAD_flowdir_land_utm.asc"
fname_outlet = regional_path + "/model/inputs/HAD_tana_dryp_station_utm.csv"
TASK: Visualize the raster datasets, use the following function to quickly create a map.
[ ]:
import matplotlib.pyplot as plt
import rasterio
# this function can be used for ploting raster maps, modify it as your conveniency.
def plot_raster_file(fname, ax=None, vmin=-20.0, vmax=20.0):
# create plot
if ax is None:
fig, ax = plt.subplots()
cmap = plt.cm.get_cmap('coolwarm_r', 12)
data = rasterio.open(fname).read(1)
im = ax.imshow(data,# origin="lower",#cmap=cmap,
#vmin=vmin, vmax=vmax,
)#extent=bounds)
ax.axis('off')
plt.colorbar(im)
return im
[ ]:
plot_raster_file(fname_surface)
[ ]:
plot_raster_file(regional_path + "/model/inputs/HAD_flowdir_D8.asc")
[ ]:
plot_raster_file(fname_flowdir)
Delineate a basin by providing a surface elevation and a the location of the basin outlet. A flow direction file can also be provided.
[ ]:
fname_out = None
[ ]:
ppbasin.get_watershed_mask(fname_surface, fname_outlet, fname_flowDir=fname_flowdir, fname_out=fname_out)
[ ]:
plot_raster_file(regional_path + "/model/inputs/HAD_DEM_utm_mm_basin.asc")
Read the file created in the previous step. This file is stored in the same directory of the surface file and it will have the same name as the surface (input file) but added “_basin” at the end (If a name of the outout filename is provided use that name in the following step).
Specify the name of the output file in the previous example.
[ ]:
fname_mask = regional_path + "/model/inputs/HAD_DEM_utm_mm_basin.asc"
[ ]:
ppbasin.get_watershed_area(fname_surface, fname_outlet, fname_flowDir=fname_flowdir, fname_mask=fname_mask)
TASK Read the csv file created in the previous state to check the catchment area
TASK: Visualize the created datasets
2. Clip basin extent by mask
We are planing to run a subbasin of the regional model, so we are going to create input files for the new model. We create this files by clipling the datasets of the regional model by the extent of the subbasin.
This following step will clip the raster files by identifiying the extend of the basin (positive values) provided as mask.
2.1. Clip flow accumulation file
[ ]:
fname = regional_path + "model/inputs/HAD_DEM_utm_mm_flowaccum.asc"
fmask = regional_path + "model/inputs/HAD_DEM_utm_mm_basin.asc"
fname_output = training_general_path + "TA_HAD_DEM_utm_mm_flowaccum.asc"
[ ]:
rrtools.clip_raster_by_mask(fname, fmask, fname_output)
2.2. Clip flow direction map
[ ]:
fname_flowdirD8 = regional_path + "model/input/HAD_flowdir_D8.asc"
fmask = regional_path + "model/input/HAD_DEM_utm_mm_basin.asc"
fname_bflowdir = basin_path + "model/input/TA_HAD_flowdir_D8_utm.asc"
[ ]:
rrtools.clip_raster_by_mask(fname_flowdirD8, fmask, fname_bflowdir)
2.3. Create model parameter datasets
TASK: Create all parameters files required for the Tana basin, make sure that all files have the same grid size. A list of files required for running the Tana basin can be found in the parameter input file.
TASK: Visualize the created datasets
3. Create a flow direction map
DRYP uses a flow direction map (optional, but recomended!) based on LandLab based notation (https://landlab.readthedocs.io/en/latest/reference/components/flow_director.html), so in order to use flow direction maps from other source that are generally in D8 format (e.g. Hydrosheds) it has to be transformed into LandLab format. The following function will help with this:
[ ]:
fflowdird8 = basin_path + "model/input/TA_HAD_flowdir_D8_utm.asc"
fname_bflowdirLandlab = basin_path + "model/input/TA_HAD_flowdir_landlab_utm.asc"
[ ]:
rrtools.create_raster_flowdirection_dryp(fflowdird8, fname_bflowdirLandlab)
TASK: Visualize the created datasets
4. Creating river network
River networks can be created form the flow accumulation map. Here, we use the flow accumalation map created in the previous step to build the river network. The river network can be created based on a minimum number of cells or the minimun area.
[ ]:
fname_flowaccum = basin_path + "model/input/TA_HAD_DEM_utm_mm_flowaccum.asc"
threshold = 10
fname_output = basin_path + "model/input/TA_HAD_river_network_utm.asc"
[ ]:
rrtools.create_raster_river_network(fname_flowaccum, threshold, fname_output, cell_area=False, fill_value=1000)
TASK: Change the threshold to see how the river network varies
TASK: Visualize the created datasets
5. Create constant head boundary condition
[ ]:
fname_wte = regional_path + "model/inputs/HAD_long_term_wte_ini.asc"
fname_outlet = regional_path + "model/inputs/HAD_tana_dryp_station_utm.csv"
fname_output = basin_path + "model/input/TA_HAD_CHB_utm.asc"
[ ]:
rrtools.create_raster_bc_at_point(fname_wte, fname_outlet, fname_output)
TASK: Modify the list of points in the csv file create a new file for head boundary conditions.
6. Clip all raster parameters files for running the model
[ ]:
import glob
files = glob.glob(regional_path+"/model/inputs/*.asc")
[ ]:
new_path = "D:/HAD_basins/Kenya/input/KE_"
fnames = []
for ifile in files:
#fnames.append(ifile.split("\\")[-1])
new_name = new_path + ifile.split("\\")[-1]
fnames.append(new_name)
print(new_name)
#new_name = new_path + fnames[0]
#process clip raster
#rrtools.clip_raster_by_mask(ifile, fmask, new_name)
7. Check that all raster datasets are correct
[ ]:
#### fname_base needs to be updated
fname_base = basin_path+"Kenya/input/KE_HAD_DEM_utm_mm.asc"
for ifname in fnames:
rrtools.check_raster_alignaments(fname_base, ifname)
[ ]: