smash.factory.detect_sink

smash.factory.detect_sink(flwdir_path, output_path=None)

Detect the sink cells in a flow direction file.

Sinks in flow direction lead to numerical issues in the routing scheme and should be removed before generating the mesh.

Parameters:

flwdir_pathstr

Path to the flow directions file. The flow direction convention is the following:

output_pathstr or None, default None

Path to the output file. If given, a raster file, in tif format, is written representing the sink cells.

Returns:

sinknumpy.ndarray

An array of shape (nrow, ncol) containing a mask of sink cells.

• 0: non-sink cell

• 1: sink cell

See also

smash.factory.generate_mesh

Automatic mesh generation.

Examples

>>> from smash.factory import load_dataset, detect_sink

Retrieve a path to a flow direction file. A pre-processed file is available in the smash package (the path is updated for each user).

>>> flwdir = load_dataset("flwdir")
flwdir

Detect the sink cells in the flow direction file.

>>> sink = detect_sink(flwdir)
>>> sink
array([[0, 0, 0, ..., 0, 0, 0],
       [0, 0, 0, ..., 0, 0, 0],
       [0, 0, 0, ..., 0, 0, 0],
       ...,
       [0, 0, 0, ..., 0, 0, 0],
       [0, 0, 0, ..., 0, 0, 0],
       [0, 0, 0, ..., 0, 0, 0]], dtype=int32)

sink is a mask of sink cells. The value 1 represents a sink cell and 0 a non-sink cell. The given flow direction file in the example does not contain any sink cell. sink is therefore a matrix of zeros.

>>> np.all(sink == 0)
np.True_

In this example, we are going to add a false sink (2 cells) to show how they can be identified.

>>> sink[10, 10:12] = 1
>>> np.all(sink == 0), np.count_nonzero(sink == 1)
(np.False_, 2)

We can retrieve the indices of the sink cells.

>>> idx = np.argwhere(sink == 1)
array([[10, 10],
       [10, 11]])

The flow direction file can be modified to remove sink cells.

>>> with rasterio.open(flwdir) as ds_in:
...     flwdir_data = ds_in.read(1)
...     flwdir_data[sink == 1] = 0
...     with rasterio.open("flwdir_wo_sink.tif", "w", **ds_in.profile) as ds_out:
...         ds_out.write(flwdir_data, 1)

Note

Setting 0 to sink cells in the flow direction file is a way to remove them. However, it might not be the best way to handle sink cells.

Finally, we can write the sink cells to a raster file by providing the output path.

>>> detect_sink(flwdir, "./flwdir_sink.tif")

The sink cells are written to the file flwdir_sink.tif and can be post-processed with a GIS software.