\cat grass
Beven et al. (1995) wrote "TOPMODEL is not a hydrological modelling package. It is rather a set of conceptual tools that can be used to reproduce the hydrological behaviour of catchments in a distributed or semi-distributed way, in particular the dynamics of surface or subsurface contributing areas."
Huidae Cho wrote r.topmodel for his Master's thesis in 2000 based on the original TOPMODEL source code (TMOD9502.FOR) and the following implementations have been developed based on his source code:
There are many combinations of input data accepted by r.topmodel. First, like many other hydrologic models, it requires a filled (depressionless) elevation map. A depressionless elevation map can be directly supplied or automatically generated by r.topmodel with the following sets of input options:
Input set 1 is usually used for the first time simulation because the user may not have various raster maps used for TOPMODEL yet. Input set 2 is good for those who already have filled elevation maps because it skips the step of creating the depressionless raster which may be the same as the one the user has given.
During the first simulation, r.topmodel generates one or two rasters and a statistics file, and these will be the same for next runs of the model. In this case, input set 3 allows the user to skip all raster map generation to save time.
Since r.topidx tries to calculate topographic index values based on the resolution of the DEM, the DEM has to have the right resolution. If the resolution is greater than the actual grid size of data, r.topidx cannot know whether cells within a grid are sinks or not.
When your elevation data is (assuming its resolution is 3 m/grid):
123 456 789
In this map, each cell composes one "data grid."
If the above data has a resolution greater than 3 m/grid (let's say 1 m/grid here), an imported raster will look like the following:
111222333 111222333 111222333 444555666 444555666 444555666 777888999 777888999 777888999
In this map, one-cell data grids have been imported as nine-cell data grids, and 9 cells compose one data grid.
Now, GRASS cannot say whether those duplicated (because of the greater resolution) 9 cells per each data grid really compose one data grid or not. Similarly, r.topidx treats those cells as sinks.
To find the correct resolution of your DEM data, run the following script (find_res.sh my_dem):
#!/bin/sh r.stats -1n "$1" | awk '{ v = $0 if(NR > 1 && v != v0){ if(min_n == 0 || n < min_n) min_n = n n = 1 v0 = v next } v0 = v n++ } END{ print min_n }'
The above script counts the minimum number of duplicated cell values. Now you want to multiply the current resolution by the printed number. For example,
> g.region -p | grep nsres | awk '{print $2}' 1 > find_res.sh my_dem 3 > g.region res=3 # 1 * 3 > d.erase > r.topidx input=my_dem output=topidx
Re: TOPIDX flow direction algorithm
r.topidx has a known issue that flat areas have topographic index values of NULL since contributing areas cannot be calculated with the D8 algorithm when area is flat. If you know any way to work around this shortcoming, please let me know.
For unit consistency over all input parameters, lengths are in meters and times are in hours except for rainfall and potential evapotranspiration time series, which are in timestep (dt).
The input file contains weather data like rainfall and potential evapotranspiration. Its format is like the following (a comment line starts with #):
# ntimesteps: Number of timesteps # dt: Time increment(timestep) # [h] # R: Rainfall # [m/dt] # Ep: Potential evapotranspiration # [m/dt] # # ntimesteps dt 66 1.0 # R Ep 0.0002 0.0 0.0003 0.0 0.0003 0.0 . . .
ntimesteps (66 in the example) is the number of R (rainfall) and Ep (potential evapotranspiration) values, so the above input file simulates 66 hours: 66 * 1 hr (ntimesteps * dt).
The format of the parameters file is as follows (descriptions of all parameters are already commented):
# Catchment name Subcatchment 1 # A: Total catchment area # [m^2] # # qs0: Initial subsurface flow per unit area # [m/h] # "The first streamflow input is assumed to represent # only the subsurface flow contribution in the watershed." # - Ref[S.C.Liaw] # lnTe: Areal average of ln(T0) = ln(Te) # [ln(m^2/h)] # m: Scaling parameter # [m] # Sr0: Initial root zone storage deficit # [m] # Srmax: Maximum root zone storage deficit # [m] # td: Unsaturated zone time delay per unit storage deficit # ( > 0.0) [h] # OR # -alpha: Effective vertical hydraulic gradient # ( <= 0.0) -10 means that alpha = 10 # vch: Main channel routing velocity # [m/h] # vr: Internal subcatchment routing velocity # [m/h] # # infex: Calculate infiltration excess if not zero (integer) # K0: Surface hydraulic conductivity # [m/h] # psi: Wetting front suction # [m] # dtheta: Water content change across the wetting front # # nch: Number of distance increments # # d: Distance from catchment # [m] # The first value should be the mainstream distance from # the subcatchment outlet to the catchment outlet. # Ad_r: Cumulative area ratio of subcatchment (0.0 to 1.0) # The first value should be zero. # # A # [m^2] 3.31697E+07 # qs0 lnTe m Sr0 Srmax td/alpha vch vr 0.000075 4. 0.0125 0.0025 0.041 60. 20000. 10000. # infex K psi dtheta 0 2. 0.1 0.1 # nch 1 # d Ad_r 8000 1.0
The original TOPMODEL (TMOD9502) requires one meta input file (topmod.run) and three data files. Topmod.run specifies three input files and the output file:
Title of the simulation inputs.dat subcat.dat params.dat topmod.out
The format of inputs.dat is almost the same as that of the input file of r.topmodel except for observed flows:
NSTEP DT R PE QOBS . . .
TMOD9502 | r.topmodel | Description |
NSTEP | ntimestep in input file | number of timesteps |
DT | dt (input file) | timestep |
R | R (input file) | rainfall |
PE | Ep (input file) | evapotranspiration |
QOBS | Qobs (option) | observed flow |
Subcat.dat contains information about the topographic index of the watershed:
NSC IMAP IOUT SUBCAT NAC AREA AC ST . . . NCH ACH D . . . MAPFILE
TMOD9502 | r.topmodel | Description |
NSC | only one subcatchment is supported per run | number of subcatchments |
IMAP | . | not implemented in TMOD9502 |
IOUT | . | level of output |
SUBCAT | Catchment name (parameters file) | subcatchment name |
NAC | nidxclass (option) | number of topographic index classes |
AREA | A (parameters file) | subcatchment area, TMOD9502: proportion of total area, r.topmodel: total area (note that r.topmodel supports only one subcatchment) |
AC | 2nd field (idxstats option) | distribution of area with topographic index, TMOD9502: sum to 1, r.topmodel: sum to total number of cells in the raster map |
ST | 1st field (idxstats option) | topographic index value |
NCH | nch (parameters file) | number of channels |
ACH | Ad_r (parameters file) | cumulative distribution of area with distance D or d from outlet |
D | d (parameters file) | distance from subcatchment outlet |
MAPFILE | . | not implemented in TMOD9502 |
Params.dat contains various information about the initial status of the model:
SUBCAT SZM T0 TD CHV RV SRMAX Q0 SR0 INFEX XK0 HF DTH
TMOD9502 | r.topmodel | Description |
SUBCAT | Catchment name (parameters file) | subcatchment name |
SZM | m (parameters file) | scaling parameter |
T0 | lnTe (parameters file) | areal average of ln(T0) |
TD | td (parameters file) | unsaturated zome time delay per unit storage deficit |
CHV | vch (parameters file) | main channel routing velocity |
RV | vr (parameters file) | internal subcatchment routing velocity |
SRMAX | Srmax (parameters file) | maximum root zone storage deficit |
Q0 | qs0 (parameters file) | initial subsurface flow per unit area |
SR0 | Sr0 (parameters file) | initial root zone storage deficit |
INFEX | infex (parameters file) | calculate infiltration excess if not zero |
XK0 | K0 (parameters file) | surface hydraulic conductivity |
HF | psi | wetting front suction |
DTH | dtheta | water content change across the wetting front |
Topmod.out contains the result of a simulation.
TMOD9502 | r.topmodel | Description |
q(it) | qt | total flow per unit area (m/timestep) |
quz | qv | vertical flux (m/timestep) |
q | qs | subsurface flow per unit area (m/timestep) |
sbar | S_mean | mean saturation deficit in the watershed (m) |
qof | qo | saturation overland flow per unit area (m/timestep) |
. | Qt | total flow (m^3/timestep) |
This example will give you a better idea of how to prepare input data. The zip file contains r.topmodel_ex/tmod9502 and r.topmodel_ex/r.topmodel directories. Basically, the two models are the same except that the former is for TMOD9502 and the latter is for r.topmodel.
Beven, K., Lamb, R., Quinn, P., Romanowicz, R., Freer, J., 1995. TOPMODEL. In: Singh, V.P. (Ed.), Computer Models of Watershed Hydrology. Water Resources Publications, pp. 627-668.