<a href="https://ofb.gouv.fr/"><img alt="Logo OFB" src="https://ofb.gouv.fr/sites/default/files/logo-ofb.png" height="50"></a> <a href="https://u-paris.fr/"><img alt="Logo UPC" src="https://u-paris.fr/wp-content/uploads/2022/03/UniversiteParisCite_logo_horizontal_couleur_RVB.jpg" height="80"></a> <a href="https://www.inrae.fr/"><img alt="Logo INRAE" src="https://seeklogo.com/images/I/inrae-logo-3BC2534099-seeklogo.com.png" height="40"></a> <a href="https://univ-cotedazur.fr/"><img alt="Logo Université Côte d'Azur" src="http://www.imbrsea.eu/sites/default/files/portfolio/UCAlogoQlarge.png" height="50"></a> <a href="https://www.ensg.eu/"><img alt="Logo ENSG" src="https://cours-fad.ensg.eu/login_ensg/img/logos/logo-ensg.png" height="50"></a> <a href="https://www.cnrs.fr/"><img alt="Logo CNRS" src="https://www.cnrs.fr/themes/custom/cnrs/logo.svg" height="70"></a>

# Diachronix plugin

## Purpose of the plugin

This plugin proposes an automatic method to segment polygon layers of valley bottom or active channel into elementary spatial units of equal length (refer here as blocks). Because valley bottom blocks are used as a reference, all the active channel blocks are consistent through time. For each of them the mean width is calculated. Thanks to this method, it is possible to track over time the active channel width evolution of specific river reaches (Figure 1), and to obtain the longitudinal pattern of active channel width along a reach.

<div align="center">
<figure>
<img src="./figures/aims.svg" alt="Aims" width="50%">
<figcaption><b>Figure 1: General purpose of the plugin.</b></figcaption>
</figure>
</div>


## Input layers

This plugin can segment valley bottom or active channel layers. The geometry of these layers has to be of *polygon* type. The coordinate reference system of these layers have to use *meters* as units. Both valley bottom and active channel layers have to be associated to a single river reach without any tributaries (Figure 2).

<div align="center">
<figure>
<img src="./figures/layers.svg" alt="Layers" width="50%">
<figcaption><b>Figure 2: Geometry of the input layers.</b></figcaption>
</figure>
</div>

As shown in Figure 2, if there are a main river and a tributary, the user has to create a valley bottom layer, or active channel layer, specifically for each river reach. The valley bottom layer can be manually digitized based on elevation raster datasets (DEM) or geological maps, or automatically extracted from available GIS tools (e.g., MRVBF from Gallant and Dowling 2003). The active channels layers can be manually digitized on ortho-images. It is better to start the digitizing by the upstream end of the valley in order to follow the flow direction.


## General method

The general procedure on which this plugin is based was designed by Frédéric Liébault and Guillaume Piton (INRAE/ETNA). The user should provide a reference layer, which in most of the cases is the valley bottom layer. In general, this layer corresponds to the spatial extent of recent alluvial deposits (active channel + floodplain + modern terraces), representing the historical channel shifting space. As previously said, this layer can be extracted from raster elevation datasets and/or geological maps. However, if it is not possible to obtain this layer, the user may just create a polygon from the fusion of all the multi-date active channel layers (from where islands are removed), to make sure that the reference layer encompasses all the available active channel extents. Once this reference layer is set, the user can digitize one, or many, active channel layers. In the case the user would like to compare active channels at different dates, one layer per date must be digitized.

Then, the reference layer is segmented in blocks of equal length. The length is provided by the user. An id and a distance from the starting point are associated to each reference layer block. After this segmentation, the active channel layers are clipped according to the reference layer blocks. The mean width is then calculated for each active channel block. Thanks to this, the block number n of the reference layer always corresponds to the block number n of the active channel layers. So the user may track over time the evolution of the width of a specific block of active channel. Figure 1 gives an overview of the method.


## The use of the plugin

### Prerequisites 

The *Diachronix* plugin requires at least QGIS 3.14. The plugin was successfully tested on Windows and on Linux. In order to be run, the *Diachronix* plugin necessitates another plugin called [Geometric attributes](https://github.com/BjornNyberg/Geometric-Attributes-Toolbox/wiki) (Nyberg et al., 2015). Moreover, the *SAGA GIS* tools should be activated. This is done automatically at the install of QGIS in most of the cases. Here, this plugin is used to calculate the centerlines of the different layers. However, the co-installation of this plugin should be automatic. Before to run the plugin the layers to be segmented should be loaded in QGIS.

### Running the plugin 

The plugin could be run either through the Extensions menu or by clicking on the corresponding icon. The main window of the plugin is divided in three steps (Figure 3).

<div align="center">
<figure>
<img src="./figures/window_parts.svg" alt="Window" width="70%">
<figcaption><b>Figure 3: The main window of the plugin.</b></figcaption>
</figure>
</div>

In step 1, the user may choose : (i) to calculate blocks of the valley bottom layer (i.e. the reference layer) only, by checking the *Valley bottom blocks* box; (ii) to calculate the blocks of the active channel layer only, by checking the *Active channel blocks* box; (iii) to calculate the blocks of both the reference layer and the active channel layer, by checking the *Valley bottom and active channel blocks* box.


### Segmentation of the reference layer

First, it is recommended to run only the segmentation of the reference layer by checking the box *Valley bottom blocks*. Once this checkbox is checked, the step 2 of the main window becomes accessible. First, in the menu *Valley bottom*, the user should select the layer corresponding to the valley bottom, or to the reference layer to be used. At this step, if we are not sure if the polygon starts from the upstream end of the valley we can click on *To detect the starting point* and *OK*. A new point layer will appear showing the starting point of the layer. However, we may note it is not always obvious to determine the starting point.

Then the user may choose a length for the blocks to be calculated at the line *Distance between two cross-sections*. By default this length is set to 100 m. The user is free to use any length but very short lengths may produce some issues later on the process. The user should then set the *Cross-sections widths* to be used. This cross-section width should be large enough to always be larger than the maximum width of the layer to be segmented.

At the line *Valley bottom blocks (output)* the user sets a path and a name for the layer containing the resulting blocks which will be created. The format of this layer may be *shapefile* or *GeoPackage*. At the line *Table of valley bottom blocks (output)* the user sets a path and a name for the table which will contain the data relative to the resulting blocks. This table can be stored in *csv*, in Excel format (*.xlsx*) or in OpenDocument format (*.ods*). This table is exactly the same as the attribute table of the produced layer. The user may choose To reverse the starting point in order to start the blocks from the upstream end of the valley or from the downstream end of the valley. In order to check the intermediate layers, the user may load the temporary layers by clicking on *To load the temporary layers*. By clicking *OK*, the blocks of the bottom valley layer will be produced. The process may take few seconds or few minutes according to the size of the layer to be segmented and the provided parameters as the distance between two cross-sections. The attribute table of the valley bottom blocks layer is shown in Figure 4.


<div align="center">
<figure>
<img src="./figures/attribute_table.png" alt="Attribute table" width="60%">
<figcaption><b>Figure 4: The attribute table of the bottom valley blocks layer.</b></figcaption>
</figure>
</div>

In this attribute table, the *fid* and *ID* fields are just two ID fields. Each block has a unique ID. The other fields are: (i) *distance*: the distance in meters between the centroid of the block and the upstream end (or the downstream end) of the valley (Figure 5); (ii) *area_m2*: the area in m<sup>2</sup> of each block; (iii) *length_m*: the length in meters of each block along its centerline; these lengths are always very close to the length set by the user for the Distance between two cross-sections; (iv) *mean_width_m*: the mean width of each block in meters; this mean width is simply the ratio of the area of each block by the length of each block.

<div align="center">
<figure>
<img src="./figures/distance_blocks.svg" alt="Distance between blocks" width="50%">
<figcaption><b>Figure 5: The length of each block, so called *distance*.</b></figcaption>
</figure>
</div>

The plugin also exports the result in a simple table format (.csv, .xlsx or .ods). This table is exactly the same as the attribute table of the exported layer but allows the user to directly work on it outside from QGIS (Figure 6).

<div align="center">
<figure>
<img src="./figures/table_ods.png" alt="Simple table" width="60%">
<figcaption><b>Figure 6: The exported table opened in LibreOffice Calc.</b></figcaption>
</figure>
</div>


### Segmentation of the active channel layer

Once the reference layer is segmented, we can use it to segment the active channel layer. To do it, in the part of the main window called *Process to be performed*, the user has to check *Active channel blocks*. By doing this, the *Active channel* part of the window becomes available. At the line *Active channel*, the user should select the layer of the active channel to be segmented. At the line *Valley bottom blocks*, the user should set the path to the reference layer to be used for the segmentation. This reference layer is the one produced at the previous step. These two layers may be in *shapefile* or in *GeoPackage*. At the line *Active channel blocks (output)*, the user should set the path and the name of the layer which will be produced after the segmentation of the active channel. At the line *Table of active channel blocks (output)*, the user should set the path and the name of the resulting table. This table could be in *.csv*, *.xlsx* or in *.ods*. Finally, it is possible to load the temporary layers by checking the line *To load the temporary layers*.

As a result, we obtain the layer of the segmented active channel. The fields of the attribute table are exactly the same as the ones of the *Bottom valley blocks*. It is the same for the exported table. Beware that the *length* parameters of these active layers might not be close to the *length* of the valley bottom layer. Indeed, a highly meandering river in a straight valley has a cumulated length quite longer than the valley axis. This length must be considered because the actual mean river width should be computed considering the cumulated channel area and its total length, including the additional distance associated with the meanders. 


### Segmentation of the reference layer and of the active channel layer in once

The user may also run the segmentation of the reference layer and of the active channel layer in once. To do this, the user should check the *Valley bottom and active channel blocks* option in the part *Process to be performed* of the main window. Once this option is checked, both the parts dedicated to the segmentation of the reference layer and to the segmentation of the active channel layer becomes available. The only change takes part in the part dedicated to the segmentation of the active channel layer. The user should not any more set the path to the segmented reference layer to be used because the process will automatically use the one produced here.


## Limitations

The plugin has some limitations listed below.


### The distance between two cross-sections

As said previously, the distance between two cross-sections is freely set by the user. However, this distance should not be too small in order to avoid some issues in the building of the blocks. The Figure 7 shows the issues when this distance is too short.

<div align="center">
<figure>
<img src="./figures/distance_too_short.png" alt="Distance too short" width="50%">
<figcaption><b>Figure 7: Effect of a distance too short (5 m in this example) on the resulting blocks.</b></figcaption>
</figure>
</div>

To avoid this issue, we have to set a larger distance between two cross-sections.


### The width of the cross-sections

The width of the cross-sections is freely set by the user. However, this width should be larger than the maximum width of the layer to be segmented. Otherwise, some blocks will be missing, resulting in blocks of non equal lengths (Figure 8).

<div align="center">
<figure>
<img src="./figures/cross-sections_too_small.svg" alt="Cross-sections too short" width="65%">
<figcaption><b>Figure 8: Effects of the cross-sections widths on the resulting blocks.</b></figcaption>
</figure>
</div>

As shown in Figure 8, the blocks 1 and 2 are not well segmented. They are two times longer than the normal blocks. This is explained by the cross-sections 1 and 2 used for the segmentation. These two cross-sections are too short compared to the maximum width of the layer to be segmented. To avoid this issue, we have to set a larger value for the width of the cross-section.


### The cross-sections and the centerlines

Cross-sections should be perpendicular to the centerlines of the layers to be segmented for a reliable computation of the channel/valley bottom width. However, this is not so obvious to implement it in an automatic method. Thus, in this plugin the cross-sections are not always perfectly perpendicular to the centerlines (Figure 9).

<div align="center">
<figure>
<img src="./figures/perpendicular.svg" alt="Cross-sections not perpendicular" width="40%">
<figcaption><b>Figure 9: The cross-section 1 is not perpendicular to the centerline.</b></figcaption>
</figure>
</div>

As shown in Figure 9, the cross-section 1 is not perpendicular to the centerline. Because of this, the area of the block 1 is slightly underestimated compared to the block 2 which is slightly overestimated. However, this limitation has relatively few impacts on the calculation of the mean widths of the blocks.


### High sinuosity

In the case of rivers with high sinuosity, some issues may appear in the meanders. In the meanders, one cross-section may cut two times the layer to be segmented. When this case happens, some blocks are missing (Figure 10).

<div align="center">
<figure>
<img src="./figures/sinuosity.svg" alt="Sinuosity" width="50%">
<figcaption><b>Figure 10: The cross-section 1 cuts two times the layer to be segmented, therefore the block a is missing, as the block b related to the cross-section 2.</b></figcaption>
</figure>
</div>

To avoid this issue, the user may shorten the width of the cross-sections or make the distance between the cross-sections larger.


## Authors, contributors and contact

The geomorphological part was carried out by [Frédéric Liébault](https://isidore.science/a/liebault_frederic) (INRAE), [Guillaume Piton](https://cv.archives-ouvertes.fr/guillaume-piton?langChosen=fr) (INRAE), [Margot Chapuis](https://cv.archives-ouvertes.fr/margot-chapuis) (Université Côte d’Azur) and [Gabriel Melun](https://www.researchgate.net/profile/Gabriel-Melun) (OFB). The implementation in Python was done by [Paul Passy](https://briques-de-geomatique.readthedocs.io/fr/latest/) (Université Paris Cité - PRODIG), [Sylvain Théry](https://art-dev.cnrs.fr/index.php/le-laboratoire/les-membres-de-l-unite/membres-permanents/212-thery-sylvain) (CNRS - ART-Dev) and students in geomatics (ENSG) : Judith Nabec, Estelle Stefanini, Damien D-Arras and Fernando Rico-Quintero


### Contact

For any comments or questions, feel free to contact the developpers: paul.passy[at]u-paris.fr, sylvain.thery[at]cnrs.fr


### References

Gallant JC, Dowling TI. 2003. A multiresolution index of valley bottom flatness for mapping depositional areas. Water Resources Research 39. DOI: [doi:10.1029/2002WR001426](https://agupubs.onlinelibrary.wiley.com/doi/full/10.1029/2002wr001426)

Nyberg B, Buckley SJ, Howell JA, Nanson RA. 2015. Geometric attribute and shape characterization of modern depositional elements: A quantitative GIS method for empirical analysis. Computers & Geosciences 82: 191-204. DOI: [https://doi.org/10.1016/j.cageo.2015.06.003](https://www.researchgate.net/profile/Bjoern-Nyberg-2/publication/279520204_Geometric_Attribute_and_Shape_Characterization_of_Modern_Depositional_Elements_A_Quantitative_GIS_Method_for_Empirical_Analysis/links/56f2aa1408ae354d162af6be/Geometric-Attribute-and-Shape-Characterization-of-Modern-Depositional-Elements-A-Quantitative-GIS-Method-for-Empirical-Analysis.pdf)

