wflow_lib Module

wflow_lib - terrain analysis and hydrological library

The goal of this module is to make a series functions to upscale maps (DEM) and to maintain as much of the information in a detailled dem when upscaling to a coarser DEM. These include:

  • river length (per cell)

  • river network location

  • elevation distribution

  • other terrain analysis

the wflow_prepare scripts use this library extensively.

$Author: schelle $ $Id: wflow_lib.py 808 2013-10-04 19:42:43Z schelle $ $Rev: 808 $

wflow_lib.Gzip(fileName, storePath=False, chunkSize=1048576)

Usage: Gzip(fileName, storePath=False, chunksize=1024*1024) Gzip the given file to the given storePath and then remove the file. A chunk size may be selected. Default is 1 megabyte Input:

fileName: file to be GZipped storePath: destination folder. Default is False, meaning the file will be zipped to its own folder chunkSize: size of chunks to write. If set too large, GZip will fail with memory problems

wflow_lib.area_percentile(inmap, area, n, order, percentile)

calculates percentile of inmap per area n is the number of points in each area, order, the sorter order of inmap per area (output of areaorder(inmap,area)) n is the output of pcr.areatotal(pcr.spatial(pcr.scalar(1.0)),area)

Input:
  • inmap

  • area map

  • n

  • order (riverorder)

  • percentile

Output:
  • percentile map

wflow_lib.area_river_burnin(ldd, dem, order, Area)

Calculates the lowest values in as DEM for each erea in an area map for river of order order

Input:
  • ldd

  • dem

  • order

  • Area map

Output:
  • dem

wflow_lib.area_riverlength_factor(ldd, Area, Clength)

ceates correction factors for riverlength for the largest streamorder in each area

Input:
  • ldd

  • Area

  • Clength (1d length of a cell (pcr.sqrt(Area))

Output:
  • distance per area

wflow_lib.areastat(Var, Area)

Calculate several statistics of Var for each unique id in Area

Input:
  • Var

  • Area

Output:
  • Standard_Deviation,Average,Max,Min

wflow_lib.checkerboard(mapin, fcc)

checkerboard create a checkerboard map with unique id’s in a fcc*fcc cells area. The resulting map can be used to derive statistics for (later) upscaling of maps (using the fcc factor)

Input:
  • map (used to determine coordinates)

  • fcc (size of the areas in cells)

Output:
  • checkerboard type map

wflow_lib.classify(inmap, lower=[0, 10, 20, 30], upper=[10, 20, 30, 40], classes=[2, 2, 3, 4])

classify a scaler maps accroding to the boundaries given in classes.

wflow_lib.configget(config, section, var, default)

Gets a string from a config file (.ini) and returns a default value if the key is not found. If the key is not found it also sets the value with the default in the config-file

Input:
  • config - python ConfigParser object

  • section - section in the file

  • var - variable (key) to get

  • default - default string

Returns:
  • string - either the value from the config file or the default value

wflow_lib.configsection(config, section)

gets the list of keys in a section

Input:
  • config

  • section

Output:
  • list of keys in the section

wflow_lib.configset(config, section, var, value, overwrite=False)

Sets a string in the in memory representation of the config object Deos NOT overwrite existing values if overwrite is set to False (default)

Input:
  • config - python ConfigParser object

  • section - section in the file

  • var - variable (key) to set

  • value - the value to set

  • overwrite (optional, default is False)

Returns:
  • nothing

wflow_lib.cutMapById(data, subcatchmap, id, x, y, FillVal)
Parameters
  • data – 2d numpy array to cut

  • subcatchmap – 2d numpy array with subcatch

  • id – id (value in the array) to cut by

  • x – array with x values

  • y – array with y values

Returns

x,y, data

wflow_lib.derive_HAND(dem, ldd, accuThreshold, rivers=None, basin=None)

Function derives Height-Above-Nearest-Drain. See http://www.sciencedirect.com/science/article/pii/S003442570800120X Input:

dem – pcraster object float32, elevation data ldd – pcraster object direction, local drain directions accuThreshold – upstream amount of cells as threshold for river

delineation

rivers=None – you can provide a rivers layer here. Pixels that are

identified as river should have a value > 0, other pixels a value of zero.

basin=None – set a boolean pcraster map where areas with True are estimated using the nearest drain in ldd distance

and areas with False by means of the nearest friction distance. Friction distance estimated using the upstream area as weight (i.e. drains with a bigger upstream area have a lower friction) the spreadzone operator is used in this case.

Output:

hand – pcraster bject float32, height, normalised to nearest stream dist – distance to nearest stream measured in cell lengths

according to D8 directions

wflow_lib.detdrainlength(ldd, xl, yl)

Determines the drainaige length (DCL) for a non square grid

Input:
  • ldd - drainage network

  • xl - length of cells in x direction

  • yl - length of cells in y direction

Output:
  • DCL

wflow_lib.detdrainwidth(ldd, xl, yl)

Determines width of drainage over DEM for a non square grid

Input:
  • ldd - drainage network

  • xl - length of cells in x direction

  • yl - length of cells in y direction

Output:
  • DCL

wflow_lib.find_outlet(ldd)

Tries to find the outlet of the largest catchment in the Ldd

Input:
  • Ldd

Output:
  • outlet map (single point in the map)

wflow_lib.getRowColPoint(in_map, xcor, ycor)

returns the row and col in a map at the point given. Works but is rather slow.

Input:
  • in_map - map to determine coordinates from

  • xcor - x coordinate

  • ycor - y coordinate

Output:
  • row, column

wflow_lib.getValAtPoint(in_map, xcor, ycor)

returns the value in a map at the point given. works but is rather slow.

Input:
  • in_map - map to determine coordinates from

  • xcor - x coordinate

  • ycor - y coordinate

Output:
  • value

wflow_lib.getcols()

returns the number of columns in the current map

Input:
Output:
  • nr of columns in the current clonemap as a scalar

wflow_lib.getgridparams()

return grid parameters in a python friendly way

Output:

[ Xul, Yul, xsize, ysize, rows, cols]

  • xul - x upper left centre

  • yul - y upper left centre

  • xsize - size of a cell in x direction

  • ysize - size of a cell in y direction

  • cols - number of columns

  • rows - number of rows

  • xlr - x lower right centre

  • ylr - y lower right centre

wflow_lib.getrows()

returns the number of rows in the current map

Input:
Output:
  • nr of rows in the current clonemap as a scalar

wflow_lib.idtoid(sourceidmap, targetidmap, valuemap)

tranfer the values from valuemap at the point id’s in sourceidmap to the areas in targetidmap.

Parameters
  • pointmap

  • areamap

  • valuemap

Returns

wflow_lib.lddcreate_save(lddname, dem, force, corevolume=1e+35, catchmentprecipitation=1e+35, corearea=1e+35, outflowdepth=1e+35)

Creates an ldd if a file does not exists or if the force flag is used

input:
  • lddname (name of the ldd to create)

  • dem (actual dem)

  • force (boolean to force recreation of the ldd)

  • outflowdepth (set to 10.0E35 normally but smaller if needed)

Output:
  • the LDD

wflow_lib.points_to_map(in_map, xcor, ycor, tolerance)

Returns a map with non zero values at the points defined in X, Y pairs. It’s goal is to replace the pcraster col2map program.

tolerance should be 0.5 to select single points Performance is not very good and scales linear with the number of points

Input:
  • in_map - map to determine coordinates from

  • xcor - x coordinate (array or single value)

  • ycor - y coordinate (array or single value)

  • tolerance - tolerance in cell units. 0.5 selects a single cell 10 would select a 10x10 block of cells

Output:
  • Map with values burned in. 1 for first point, 2 for second and so on

wflow_lib.pt_flow_in_river(ldd, river)

Returns all points (True) that flow into the mak river (boolean map with river set to True)

Parameters
  • ldd – Drainage network

  • river – Map of river (True River, False non-river)

Return ifmap

map with infrlo points into the river (True)

Return ctach

catchment of each of the inflow points

wflow_lib.readMap(fileName, fileFormat)

Read geographical file into memory

Parameters
  • fileName

  • fileFormat

Return x, y, data, FillVal

wflow_lib.riverlength(ldd, order)

Determines the length of a river using the ldd. only determined for order and higher.

Input:
  • ldd, order (streamorder)

Returns:
  • totallength,lengthpercell, streamorder

wflow_lib.sCurve(X, a=0.0, b=1.0, c=1.0)

sCurve function:

Input:
  • X input map

  • C determines the steepness or “stepwiseness” of the curve. The higher C the sharper the function. A negative C reverses the function.

  • b determines the amplitude of the curve

  • a determines the centre level (default = 0)

Output:
  • result

wflow_lib.sCurveSlope(X, a=0.0, b=1.0, c=1.0)

First derivative of the sCurve defined by a,b,c at point X

Input:
  • X - value to calculate for

  • a

  • b

  • c

Output:
  • first derivative (slope) of the curve at point X

wflow_lib.snaptomap(points, mmap)

Snap the points in _points_ to nearest non missing values in _mmap_. Can be used to move gauge locations to the nearest rivers.

Input:
  • points - map with points to move

  • mmap - map with points to move to

Return:
  • map with shifted points

wflow_lib.subcatch(ldd, outlet)

Determines a subcatchment map using LDD and outlet(s). In the resulting subcatchment map the i’s of the catchment are determiend by the id’s of the outlets.

Input:
  • ldd

  • Outlet - maps with points for each outlet.

Output:
  • map of subcatchments

wflow_lib.subcatch_order_a(ldd, oorder)

Determines subcatchments using the catchment order

This version uses the last cell BELOW order to derive the catchments. In general you want the _b version

Input:
  • ldd

  • order - order to use

Output:
  • map with catchment for the given streamorder

wflow_lib.subcatch_order_b(ldd, oorder, sizelimit=0, fill=False, fillcomplete=False, stoporder=0)

Determines subcatchments using the catchment order

This version tries to keep the number op upstream/downstream catchment the small by first dederivingatchment connected to the major river(the order) given, and fill up from there.

Input:
  • ldd

  • oorder - order to use

  • sizelimit - smallest catchments to include, default is all (sizelimit=0) in number of cells

  • if fill is set to True the higer order catchment are filled also

  • if fillcomplete is set to True the whole ldd is filled with catchments.

:returns sc, dif, nldd; Subcatchment, Points, subcatchldd

wflow_lib.subcatch_stream(ldd, threshold, min_strahler=- 999, max_strahler=999, assign_edge=False, assign_existing=False, up_area=None)

(From Deltares Hydrotools)

Derive catchments based upon strahler threshold Input:

ldd – pcraster object direction, local drain directions threshold – integer, strahler threshold, subcatchments ge threshold

are derived

min_strahler – integer, minimum strahler threshold of river catchments

to return

max_strahler – integer, maximum strahler threshold of river catchments

to return

assign_unique=False – if set to True, unassigned connected areas at

the edges of the domain are assigned a unique id as well. If set to False, edges are not assigned

assign_existing=False == if set to True, unassigned edges are assigned

to existing basins with an upstream weighting. If set to False, edges are assigned to unique IDs, or not assigned

output:

stream_ge – pcraster object, streams of strahler order ge threshold subcatch – pcraster object, subcatchments of strahler order ge threshold

wflow_lib.sum_list_cover(list_of_maps, covermap)

Sums a list of pcrastermap using cover to fill in missing values

Parameters
  • list_of_maps – list of maps to sum

  • covermap – maps/ value to use fro cover

Returns

sum of list of maps (single map)

wflow_lib.upscale_riverlength(ldd, order, factor)

Upscales the riverlength using ‘factor’ The resulting maps can be resampled (e.g. using resample.exe) by factor and should include the accurate length as determined with the original higher resolution maps. This function is depricated, use are_riverlength instead as this version is very slow for large maps

Input:
  • ldd

  • minimum streamorder to include

Output:
  • distance per factor cells

wflow_lib.writeMap(fileName, fileFormat, x, y, data, FillVal)

Write geographical data into file

wflow_lib.zipFiles(fileList, fileTarget)

Usage: zipFiles(fileList, fileTarget) zip the given list of files to the given target file Input:

fileList: list of files to be zipped fileTarget: target zip-file