.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "example_gallery/cell_centric_operations/downslope_path.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code. .. rst-class:: sphx-glr-example-title .. _sphx_glr_example_gallery_cell_centric_operations_downslope_path.py: Downslope path ============== Use a DEM to determine the downslope path from a point on the DEM to the edge. DEM data source: https://www.eea.europa.eu/data-and-maps/data/copernicus-land-monitoring-service-eu-dem Introduction ------------ In this example the downslope path is found from a given point on a DEM. This can be conceptualized as the path water would flow donw from a given source. This algorithm does not stop at a local minimum, but continues when the edge of the DEM is reached. The path can therefore be best conceptualized as a stram, possibly with small lakes, rather than a single droplet. Method ------ This is done by iteratively finding the neighbour with the lowest elevation. In this example, the flow is as follows: 1. A cell is supplied to the "downslope_path" function 2. The neighbouring cells are obtained 3. These cells are added to a list containing all neighbours along the path so far 4. The neighbour along the path so far with the lowest elevation is added to the path, at which point we start again at 1. The loop stops when NaNs are encountered, which happens when the border of the DEM crop is reached. This is of course assuming there are no nodata_values in the DEM. First let's import the dependencies and read in the DEM. .. GENERATED FROM PYTHON SOURCE LINES 34-45 .. code-block:: Python import matplotlib.pyplot as plt import numpy from gridkit import read_raster dem = read_raster( "../../tests/data/alps_dem.tiff", bounds=(29200, 167700, 29800, 168100) ) .. GENERATED FROM PYTHON SOURCE LINES 47-52 Next, set up the recursive function and the required empty lists. .. Warning :: Python's recursion limit can be reached with this approach for paths longer than shown here. While you can raise this limit using `sys.setrecursionlimit`, you may want to reconsider the approach and use a while-loop instead of recursion. .. GENERATED FROM PYTHON SOURCE LINES 53-81 .. code-block:: Python visited_cells = [] all_neighbours = [] neighbour_values = [] def downslope_path(cell_id): visited_cells.append(cell_id) # add new neigbours and their values to their respective lists neighbours = dem.neighbours(cell_id, connect_corners=True) for cell, value in zip(neighbours, dem.value(neighbours)): cell = tuple(cell.index) if not cell in visited_cells and not cell in all_neighbours: all_neighbours.append(cell) neighbour_values.append(value) # stop when we encounter nodata values, assuming these only occur byound the bounds of the DEM crop if dem.nodata_value in neighbour_values: return # determine the next cell and remove it from the neighbour-lists min_id = numpy.argmin(neighbour_values) next_cell = all_neighbours[min_id] all_neighbours.pop(min_id) neighbour_values.pop(min_id) return downslope_path(next_cell) .. GENERATED FROM PYTHON SOURCE LINES 82-83 Determine the starting cell and call the function .. GENERATED FROM PYTHON SOURCE LINES 83-87 .. code-block:: Python start_cell_id = dem.cell_at_point((29330, 167848)) downslope_path(tuple(start_cell_id.index)) .. GENERATED FROM PYTHON SOURCE LINES 88-93 Results ------------- Let's plot the path on the DEM For each cell along the path we plot the polygon as a square. This looks better than points at the centroids. .. GENERATED FROM PYTHON SOURCE LINES 94-114 .. code-block:: Python im = plt.imshow(dem, extent=dem.mpl_extent) plt.colorbar(fraction=0.032, pad=0.01, label="Elevation [m]") plt.xlabel("longitude") plt.ylabel("latitude") plt.title("Downslope path from a chosen starting point", fontsize=10) path_label = "Downslope path" for poly in dem.to_shapely(visited_cells): x, y = poly.exterior.xy plt.fill(x, y, alpha=1.0, color="red", label=path_label) path_label = None # Only show the label in the legend once plt.scatter( *dem.centroid(start_cell_id), marker="*", label="Starting point", color="purple" ) plt.legend() plt.show() .. image-sg:: /example_gallery/cell_centric_operations/images/sphx_glr_downslope_path_001.png :alt: Downslope path from a chosen starting point :srcset: /example_gallery/cell_centric_operations/images/sphx_glr_downslope_path_001.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 115-118 While the path is merely visualized in this example, the obtained cells can of course be used in different ways. They can for example be used to obtain values form a different dataset along the path, or intersected with a polygon to determine if the path crosses a particular zone of interest. .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 1.126 seconds) .. _sphx_glr_download_example_gallery_cell_centric_operations_downslope_path.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: downslope_path.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: downslope_path.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: downslope_path.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_