The results revealed a significant increase with time in the number of events that produce lower bed incision and bed aggradation. The number of minor events also increased. However, the frequency of bankfull events, which are responsible for dominant scouring processes, tended to decrease.
4.1. Retrospective Simulation of Peak Flows and Potential Sediment Inputs Using GeoWEPP
The GeoWEPP model has made it possible to delimit the catchment areas of the two study RCRs and to trace their respective drainage networks. In addition, it has provided information regarding runoff, soil loss, deposition, and sediment production from slopes and stream stretches. Generally, the results of the evaluation are mapped as a measurement relative to a tolerable soil loss or standard value (T), where T = 11.2 t ha-1 year-1. The non-tolerable rates, located above the T threshold, are mainly concentrated in the upper sub-catchment and headwater areas, where steep slopes (>25%) devoid of vegetation concur with the presence of quite fragile metamorphic materials (slates, phyllites, and shales). When comparing the sediment yield maps obtained at the sub-watershed level for five-year periods (
Figure 4), two clear phases can be distinguished: a first stage (1996-2012) with high specific soil degradation and intense delivery of sediments, in which the hillslopes with CSY and ASY above 150 and 20 t, respectively, abound; and another (2012-2022), characterized by more moderate and low sediment production rates, coinciding with lower runoff coefficients (
Figure 4).
Figure 4.
Spatial distribution of the cumulative (CSY) and average (ASY) sediment yields, estimated using GeoWEPP for retroactive five-year periods in the middle and upper sub-catchments of the Azohía Rambla.
Figure 4.
Spatial distribution of the cumulative (CSY) and average (ASY) sediment yields, estimated using GeoWEPP for retroactive five-year periods in the middle and upper sub-catchments of the Azohía Rambla.
In
Figure 5 simulation values of maximum runoff and potential sediment inputs estimated by the watershed and flow paths methods of WEPP are graphically displayed for each RCR at the event scale, while
Table 6 shows the mean values predicted for each PECP. The temporal distribution of both the peak discharges and the delivery of sediments for transport follows two clear patterns: 1) one, represented by the period 2010-2022, which is characterized by a high occurrence of class 2 events (0.46 times/year), corresponding to bankfull discharges, and a low presence of class 1, 3, and 4 events (0.08, 0.17, and 0.25 times/year, respectively); and 2) a second pattern, observed in the stage 1996-2009, with higher frequencies of all classes of events: 0.21, 0.50, 0.50, and 0.57 times/year, implying shorter return periods.
The SY values show similar patterns at the event level, in both UPR and MDR. The greatest differences between the two sections occur for class 1 events, during large flash floods, in which the MDR is exposed to potential inputs of 9,000 to 27,000 t of sediments, compared to 4,700-11,800 t generated by erosion in the sub-watersheds draining to UPR. This represents an increase of 90 to 137% in the amount of potentially transportable sediments along the MDR with respect to the UPR. This increase is reduced as the category of events decreases: between 80 and 115% for EC2, and below 80% for EC3 and EC4. It should be noted that, although most of the sediments are transported in suspension after superficial soil washing processes in the sub-watershed areas, one of the major morphodynamic drivers in this type of gravel-bed stream is, instead, the bedload transport [
44]. For both RCRs, the relationships between bedload transport and water discharge are dominated by large particles and sudden changes in excess energy during isolated peak discharges in time. This behavior is consistent with the results found by numerous authors [
4,
5,
45], who classified the bankfull or sub-bankfull discharge as effective discharge, in our case characterized by a recurrence interval of one to two years.
Figure 5.
Peak discharges and potential average sediment outputs estimated with GeoWEPP for the flow events that occurred during the period 1996-2022.
Figure 5.
Peak discharges and potential average sediment outputs estimated with GeoWEPP for the flow events that occurred during the period 1996-2022.
In relation to the PECPs, a clear contrast exists between the 1996-2009 and 2010-2022 stages. However, in this case, the patterns of hydrological and morphosedimentary behavior in both RCRs are defined by different associations and frequencies of ECPs. The most recent period (2010-2022) is characterized by a high recurrence of D patterns made up of a single class 2 event and the confirmation of pattern A, also endowed with a single event, this time of the flash flood type and overtopping stage. For its part, the period 1996-2009 shows a different sequential pattern dominated by pattern B, which groups class 2, 3, and 4 events (
Table 6).
Table 6.
Simulation with GeoWEPP of the potential sediment inputs entering the UPR and MDR reaches for the periods established, according to event class patterns (ECP).
Table 6.
Simulation with GeoWEPP of the potential sediment inputs entering the UPR and MDR reaches for the periods established, according to event class patterns (ECP).
PECP |
EC |
ECP |
UPR |
MDR |
ASIP (t) |
CSIP (t) |
ASIP (t) |
CSIP (t) |
July 2020 – Feb. 2023 |
2 |
D |
2208.1 |
2208.1 |
4222.0 |
4222.0 |
Sep. 2019 - July 2020 |
3-4-2 |
B |
1907.6 |
3815.1 |
2708.2 |
9694.8 |
Sep. 2018 - Sep.2019 |
1 |
A |
4712.4 |
4712.4 |
9694.8 |
9694.8 |
Oct. 2013 – Sep. 2018 |
4-4-3-2 |
B |
1181.2 |
4724.7 |
2112.8 |
2112.8 |
Sep. 2012 – Oct. 2013 |
2 |
D |
2760.7 |
2760.7 |
5000.6 |
5000.6 |
Aug. 2012 – Sep. 2012 |
2 |
D |
2200.1 |
2200.1 |
4364.0 |
4364.0 |
June 2010 – Aug. 2012 |
2 |
D |
1977.0 |
1977.0 |
3889.2 |
3889.2 |
March 2009 – June 2010 |
4-1-2-3 |
C |
3343.8 |
13375.4 |
7144.4 |
28577.6 |
Oct. 2007 – March 2009 |
2-3 |
B |
2159.9 |
4319.8 |
3963.9 |
7927.9 |
Nov. 2005 – Oct. 2007 |
2-1-3-4 |
C |
3411.1 |
10233.3 |
6082.3 |
18246.8 |
Nov. 2003 – Nov. 2005 |
3-4-2 |
B |
1670.5 |
5011.6 |
3349.8 |
10049.5 |
March 2002 – Nov. 2003 |
3-2-4 |
B |
1706.8 |
5120.5 |
3243.5 |
9730.5 |
Dec. 1998 - March 2002 |
4-4-1-3-3 |
E |
3262.1 |
16310.5 |
6990.3 |
34951.6 |
Oct. 1997 – Dec. 1998 |
2-4-3 |
B |
1500.7 |
4502.1 |
2826.0 |
8478.0 |
Oct. 1996 – Oct. 1997 |
2 |
D |
2612.0 |
2612.0 |
4711.0 |
4711.0 |
In terms of mean values, the potential inputs of sediments per event do not differ appreciably between the patterns described in UPR, especially between B and D, but in MDR more distinguishing features are observed, particularly when comparing the PECPs with patterns that include class 1 events with those that do not. More pronounced are the differences between the CSIP values of the PECPs represented by patterns A, C, and E and those of type B and D. Without a doubt, the largest potential cumulative sediment inputs registered at the entries into both sections correspond to the PECPs with patterns C (November 2005 – October 2007; March 2009 – June 2010) and E (December 1998 – March 2002), which include overflow and bankfull peak discharges. These caused important overall morphological adjustments that affected both the channel geometry and the bedform units.
4.2. Retrospective Changes in Bed Elevation and Sediment Budgets Derived from SfM-MVS Generated HRDTMs and the Resulting DoDs
The results for other time periods derived from experimental models simulating channel morphology have been found to be weakly dependent on the grid resolution [
46,
47,
48]. We used HRDTMs here because the poor resolution of the other grid models has proven in many cases to be insufficient to quantify net sediment fluxes or unit volumes of bed erosion and deposition, or to represent well the different morphosedimentary features and their changes over time. Lower resolution DTMs were ruled out due to the low geometric and textural definition they offer in unit bedform identification. This is the case of the resulting NPAO DTMs and DoDs (
Figure 6), whose loss of precision with respect to the SfM and TLS products would significantly worsen the results of a retrospective morphodynamic simulation. Nicholas et al. (2013) examined the model sensitivity to grid resolution, morphodynamic scaling, inlet boundary conditions, and parameterization of bed roughness and sediment transport, and reached similar conclusions for their simulations.
The main influence of grid resolution is on bar size, because grid resolution sets a limit on the smallest features that can be considered. Consequently, the finest grid (pixel size 0.02 m) used in our retrospective simulation allows the recognition of the smallest bars and detailed monitoring of their structural and geometric changes. Overall, the modeled channel width distributions from NPAO DTMs are consistent with the field data, although there is a trend for secondary channels and the smallest alluvial bars to be underrepresented in the resulting channel morphology relative to those monitored using the HRDoDs.
Clearly, the controls on unit bar frequency and runs in this EGBS are complex and likely include very different morphological responses of the channel to intrinsic fluvial processes, depending on the event class and the class association patterns. The detection of these morphological adjustments requires a combination of high-resolution 3D models and ground-based monitoring on the flow regime and sediment transport mechanics. In principal, unlike channels with a high width-to-depth ratio and extensive development of braided bars, where the results of two-dimensional and three-dimensional model simulations are comparable [
49], in our case the HRDoDs do adequately resolve near-bed flow, since they were implemented with very fine spatial resolutions [
50].
Figure 6.
DoDs resulting from the difference between the DEMs derived from the NPAO 2009 and 2016 LiDAR point clouds (0.5 m per pixel) and those generated with SfM-MVS in 2018 and 2023 after being converted to the same resolution. Upper and Middle RCRs.
Figure 6.
DoDs resulting from the difference between the DEMs derived from the NPAO 2009 and 2016 LiDAR point clouds (0.5 m per pixel) and those generated with SfM-MVS in 2018 and 2023 after being converted to the same resolution. Upper and Middle RCRs.
The morphological sediment budget data obtained from DoD with grid resolution of 0.5 m/pixel, using NPAO DEMs and SfM datasets, are much less realistic than those provided solely by SfM-MVS photogrammetry during the monitoring period (2018-2022). For this period the % error attributed to ANTD, UVSL, and UVSR in both RCRs varies between 0.112 and 0.192 for the resolution of 0.5 m per pixel of the NPAO DoDs and between 0.047 and 0.098 (
Table 7) when using the HRDoDs of SfM (
Table 9).
Table 7.
Statistical descriptors relating to the morphological sediment budgets calculated for the reference channel reaches (RCRs) in the UPR and MDR, for the monitoring (1) and simulation (2) periods (results from DoD generated using NPAO DEMs and SfM datasets).
Table 7.
Statistical descriptors relating to the morphological sediment budgets calculated for the reference channel reaches (RCRs) in the UPR and MDR, for the monitoring (1) and simulation (2) periods (results from DoD generated using NPAO DEMs and SfM datasets).
|
TAI |
TNVD |
ANTD |
PI |
TASL |
TASR |
UVSL |
UVSR |
SD * |
Period |
RCR |
m2
|
m3
|
Error (p.u.) |
m |
Error (p.u.) |
(p.u) |
m2
|
m2
|
m3 m-2
|
Error (p.u.) |
m3 m-2
|
Error (p.u.) |
m |
2016- 2022 (1)
|
UPR |
2407 |
440.1 |
0.112 |
0.183 |
0.112 |
0.442 |
155.3 |
2251.4 |
0.185 |
0.169 |
0.208 |
0.192 |
0.188 |
MDR |
3956 |
692.3 |
0.175 |
0.175 |
0.175 |
0.459 |
364.0 |
3591.9 |
0.085 |
0.115 |
0.201 |
0.189 |
0.106 |
2009 -2016(2)
|
UPR |
2435 |
-366.7 |
-0.089 |
-0.151 |
-0.089 |
-0.448 |
2261.2 |
173.7 |
0.172 |
0.196 |
0.122 |
0.174 |
0.156 |
MDR |
4040 |
-501.4 |
-0.106 |
-0.124 |
-0.106 |
-0.463 |
3743.6 |
296.8 |
0.139 |
0.206 |
0.067 |
0.156 |
0.116 |
In addition, the mean, total, and unit values of bed morphological change indicators that were calculated from the NPAO DoDs are quite different from those provided by the SfM HRDoDs in all of the reaches, for both the monitoring and simulation periods. In
Table 8 it can be verified that the dissimilarity coefficients between the results extracted from both types of DoDs are mostly above 0.15 or -0.15, which means that the values deduced from the lower resolution model overestimate or underestimate by more than 15%, respectively, those monitored with UAVs. Only the TAI and TASR values in both RCRs, PI in UPR, and UVSL in MDR present dissimilarity coefficients ≤ 0.15.
Table 8.
Dissimilarity coefficients between the DoDs derived from the PNOA DEMs and those obtained using SfM-MVS photogrammetry for the monitoring (1) and simulation (2) periods in UPR and MDR.
Table 8.
Dissimilarity coefficients between the DoDs derived from the PNOA DEMs and those obtained using SfM-MVS photogrammetry for the monitoring (1) and simulation (2) periods in UPR and MDR.
Period |
RCR |
TAI |
TNVD |
ANTD |
PI |
TASL |
TASR |
UVSL |
UVSR |
SD* |
2016 – 2022 (1)
|
UPR |
0.09 |
-0.23 |
-0.30 |
0.01 |
0.47 |
0.04 |
-0.64 |
-0.26 |
-0.43 |
MDR |
0.15 |
-0.33 |
-0.43 |
-0.18 |
0.60 |
0.04 |
0.00 |
-0.28 |
0.17 |
2009 - 2016 (2)
|
UPR |
0.04 |
-0.16 |
-0.22 |
0.09 |
-0.64 |
0.69 |
-0.12 |
0.16 |
-0.37 |
MDR |
0.04 |
-0.88 |
-0.32 |
-0.51 |
-0.53 |
0.79 |
0.11 |
0.66 |
0.34 |
Consequently, for the case study, where the small bed forms within a narrow ephemeral gravel-bed channel play an important role in its general morphological evolution, it is necessary to use HRDoDs to better quantify the changes produced at the event scale and PECPs. Bed adjustments at this temporal level are usually very varied and complex; in some cases, scour and downcutting phenomena predominate, while in others deposition and entrainment processes and superficial washing are prevalent. Using HRDoDs, therefore, it is possible to more precisely detect the current changes produced in this type of stream, including minor adjustments, and to obtain higher quality simulations of short-term retrospective morphodynamic changes. A good indicator of these adjustments is the sediment budget in the reaches with the highest bedload. In particular, in the RCRs analyzed here, the bedload experienced important spatial variations at the event scale, significantly affecting the sediment sources (areas of erosion) and sinks (areas of deposition) in the period 2018-2020[
4]. Once the monitoring period has been extended to 2022, it was possible to establish grouping patterns of event classes that allowed us to go back in time looking for their occurrence. Fitted HRDoDs were applied to each of them according to their degree of similarity with the starting ECPs in order to simulate their possible changes in bed elevation and sediment budgets during the period 1996-2018, the results of which are shown in
Figure 7 and
Table 9.
Figure 7.
Retrospectively simulated HRDoDs for five-year sub-periods from 1997 to 2017, combining ECPs with UAV-SfM-derived DTMs over the monitoring stage (2018-2022).
Figure 7.
Retrospectively simulated HRDoDs for five-year sub-periods from 1997 to 2017, combining ECPs with UAV-SfM-derived DTMs over the monitoring stage (2018-2022).
As already indicated in section 3.5, a total of five events were monitored and classified into four different classes: a flash flood –Class 1- (19-20 April, 2019), two bankfull events –Class 2- (24 March, 2020 and 4-5 April, 2022), a sub-bankfull flow –Class 3- (12 September, 2019), and a minor event, categorized as sub-half-bankfull –Class 4- (20 January, 2020), all having different geomorphic impacts. In the September 2018–September 2019 stage, lateral erosion from steep alluvial banks, active low bars, partially destroyed coarse bar heads, and finer-grained bar tails along the upstream reaches provided a large bedload in the downstream direction. As a result, the greatest deposition thicknesses were recorded in the flanks of the longitudinal and medial alluvial bars, in both RCRs [
4]. Previous field observations suggest that during low-water stages, as in the event of November 2018, when the top of the bar emerged, vertical accretion of these bars ceases and new secondary channels form, causing small island bars to migrate. The flash flood of 19–20 April, 2019 resumed the aggradation process, with very widespread increases in bed height. The PECP defined by the pattern of this class 1 event, when it occurs in isolation or is preceded by events of insignificant impact, was classified as a type A PECP.
Table 9.
Statistical descriptors relating to the morphological sediment budgets calculated for the reference channel reaches (RCRs) during the period monitored with SfM data (September 2018 to February 2023) and the periods simulated with ECP (October 1996 to September 2018).
Table 9.
Statistical descriptors relating to the morphological sediment budgets calculated for the reference channel reaches (RCRs) during the period monitored with SfM data (September 2018 to February 2023) and the periods simulated with ECP (October 1996 to September 2018).
|
TAI |
TNVD |
ANTD |
PI |
TASL |
TASR |
UVSL |
UVSR |
*SD |
PECP |
RCR |
m2
|
m3
|
Error (u.p.) |
m |
Error (u.p.) |
p.u. |
m2
|
m2
|
m3 m-2
|
Error (p.u.) |
m3 m-2
|
Error (p.u.) |
m |
JUL20- FEB23 |
UPR |
2730 |
395.2 |
0.063 |
0.145 |
0.063 |
0.477 |
81.2 |
2649.3 |
0.120 |
0.057 |
0.153 |
0.079 |
0.136 |
MDR |
4532 |
118.7 |
0.054 |
0.026 |
0.054 |
0.069 |
1266.0 |
3266.2 |
0.294 |
0.086 |
0.150 |
0.113 |
0.146 |
SEP19- JUL20 |
UPR |
2976 |
-575.3 |
-0.048 |
-0.193 |
-0.048 |
-0.420 |
2744.1 |
232.2 |
0.229 |
0.043 |
0.237 |
0.040 |
0.162 |
MDR |
4707 |
-742.2 |
-0.055 |
-0.158 |
-0.055 |
-0.319 |
4043.8 |
663.4 |
0.235 |
0.042 |
0.317 |
0.030 |
0.243 |
SEP18- SEP19 |
UPR |
2763 |
613.2 |
0.044 |
0.222 |
0.044 |
0.486 |
67.3 |
2695.1 |
0.128 |
0.070 |
0.231 |
0.043 |
0.118 |
MDR |
4885 |
1013.1 |
0.046 |
0.207 |
0.046 |
0.486 |
168.1 |
4717.0 |
0.086 |
0.102 |
0.218 |
0.046 |
0.106 |
OCT13- SEP18 |
UPR |
2551 |
-680.5 |
-0.112 |
-0.267 |
-0.112 |
-0.494 |
2516.7 |
34.7 |
0.272 |
0.198 |
0.127 |
0.109 |
0.105 |
MDR |
4490 |
-716.1 |
-0.106 |
-0.159 |
-0.106 |
-0.306 |
3873.6 |
616.9 |
0.243 |
0.206 |
0.368 |
0.185 |
0.186 |
SEP12- OCT13 |
UPR |
2543 |
443.8 |
0.089 |
0.174 |
0.089 |
0.485 |
51.9 |
2491.4 |
0.129 |
0.095 |
0.181 |
0.099 |
0.077 |
MDR |
4155 |
336.3 |
0.093 |
0.081 |
0.093 |
0.272 |
1032.9 |
3122.0 |
0.136 |
0.112 |
0.153 |
0.123 |
0.164 |
AUG12- SEP12 |
UPR |
2543 |
395.3 |
0.146 |
0.155 |
0.146 |
0.485 |
51.9 |
2491.4 |
0.115 |
0.122 |
0.161 |
0.148 |
0.069 |
MDR |
4155 |
298.0 |
0.157 |
0.072 |
0.157 |
0.272 |
1032.9 |
3122.0 |
0.121 |
0.158 |
0.135 |
0.141 |
0.146 |
JUN10 – AUG12 |
UPR |
2543 |
356.8 |
0.167 |
0.140 |
0.167 |
0.485 |
51.9 |
2491.4 |
0.104 |
0.186 |
0.145 |
0.109 |
0.062 |
MDR |
4155 |
270.7 |
0.189 |
0.065 |
0.189 |
0.272 |
1032.9 |
3122.0 |
0.110 |
0.199 |
0.123 |
0.098 |
0.132 |
MAR09 JUN10 |
UPR |
2542 |
-755.9 |
-0.129 |
-0.297 |
-0.129 |
-0.494 |
2508.5 |
33.9 |
0.303 |
0.202 |
0.146 |
0.133 |
0.118 |
MDR |
4155 |
-894.3 |
-0.187 |
-0.215 |
-0.187 |
-0.434 |
3746.1 |
408.5 |
0.257 |
0.274 |
0.166 |
0.145 |
0.176 |
OCT07- MAR09 |
UPR |
2542 |
-78.5 |
-0.196 |
-0.031 |
-0.196 |
-0.120 |
1438.2 |
1104.2 |
0.141 |
0.143 |
0.112 |
0.089 |
0.177 |
MDR |
4155 |
112.4 |
-0.213 |
0.027 |
-0.213 |
0.074 |
1900.3 |
2254.4 |
0.169 |
0.153 |
0.193 |
0.149 |
0.252 |
NOV05- OCT07 |
UPR |
2542 |
-802.9 |
-0.229 |
-0.316 |
-0.229 |
-0.494 |
2508.5 |
33.9 |
0.322 |
0.224 |
0.156 |
0.144 |
0.125 |
MDR |
4155 |
-983.7 |
-0.258 |
-0.237 |
-0.258 |
-0.434 |
3746.1 |
408.5 |
0.282 |
0.239 |
0.182 |
0.175 |
0.194 |
NOV03- NOV05 |
UPR |
2542 |
-64.2 |
-0.216 |
-0.025 |
-0.216 |
-0.120 |
1438.2 |
1104.2 |
0.115 |
0.106 |
0.092 |
0.107 |
0.145 |
MDR |
4155 |
85.0 |
-0.224 |
0.020 |
-0.224 |
0.074 |
1900.3 |
2254.4 |
0.128 |
0.125 |
0.146 |
0.127 |
0.191 |
MAR02- NOV03 |
UPR |
2533 |
-496.9 |
-0.285 |
-0.196 |
-0.285 |
-0.407 |
2128.1 |
405.0 |
0.260 |
0.259 |
0.140 |
0.121 |
0.219 |
MDR |
4153 |
-771.8 |
-0.297 |
-0.186 |
-0.297 |
-0.455 |
3845.1 |
308.1 |
0.211 |
0.189 |
0.125 |
0.095 |
0.160 |
DEC98- MAR02 |
UPR |
2542 |
-847.4 |
-0.219 |
-0.333 |
-0.219 |
-0.494 |
2508.5 |
33.9 |
0.340 |
0.287 |
0.164 |
0.133 |
0.132 |
MDR |
4155 |
-987.5 |
-0.246 |
-0.238 |
-0.246 |
-0.434 |
3746.1 |
408.5 |
0.284 |
0.256 |
0.183 |
0.148 |
0.194 |
OCT97- DEC98 |
UPR |
2533 |
-754.6 |
-0.253 |
-0.298 |
-0.253 |
-0.494 |
2502.2 |
30.9 |
0.303 |
0.288 |
0.142 |
0.185 |
0.116 |
MDR |
4153 |
-879.2 |
-0.267 |
-0.212 |
-0.267 |
-0.434 |
3745.1 |
408.0 |
0.253 |
0.249 |
0.163 |
0.149 |
0.173 |
OCT96- OCT97 |
UPR |
2533 |
421.5 |
0.267 |
0.166 |
0.267 |
0.486 |
50.0 |
2483.1 |
0.121 |
0.166 |
0.172 |
0.179 |
0.073 |
MDR |
4153 |
319.2 |
0.288 |
0.077 |
0.288 |
0.272 |
1032.4 |
3120.7 |
0.129 |
0.172 |
0.145 |
0.195 |
0.156 |
The statistical descriptors of the morphological sediment budgets in
Table 9 show significant cumulative changes in ground surface elevation after this monitored A stage. An average net thickness difference of around 22 and 21 cm was found in the upper and middle reaches, respectively. The sedimentary balance was positive in both RCRs, verifying a clear dominance of deposition over erosion. The proximity of both reaches to abundant sources of coarse sediment together with a locally strong connectivity between the channel bed and active alluvial banks and almost similar bed slopes explain the slight differences in the average sediment budget between UPR and MDR [
4]. However, the high average deposition unit rate observed in the middle reach (UVSL = 0.22 m3 m−2) was not accompanied by equally significant unit volumes of erosion upstream (UVSR in UPR = 0.13 m3 m−2). In addition, the net unit bed erosion rate found in MDR (UVSL < 0.09 m3 m−2) was lower than for UPR (
Table 9). Thus, the bed deposition in a downstream direction could only progress due to bank erosion and bedload mobility along intermediate sections between the two RCRs, where bank breaking phenomena and gravel removal phenomena are quite recurrent. It should be taken into account that the average unit volumes of surface lowering or bed growth are referred to the changing area, so that the highest values of UVSL or UVSR do not always correspond to higher bed material transport rates and larger volumetric change. In fact, if we consider the PECP September 2018 to September 2019 as an example, in which a major class 1 event occurred, UVSR was higher in UPR (0.231 m3 m-2) than in MDR (0.218 m3 m-2). However, since the affected area was much larger in the latter reach (4717 m2 compared to 2695 m2 in UPR) (
Table 9), the total volume of vertical bed accretion in MDR (1028 m3) almost doubled that registered upstream along the upper reach (622 m3). When retrospectively analyzing the entire period, two clear phases were distinguished: 1) an older one (1997-2012), characterized by dominant incision processes in both sections that accounted for a total excavated volume of 2,626 and 3,434 m3 in UPR and MDR, respectively; and 2) well-defined channel a more recent one (2012-2022), in which the sedimentation of coarse material and the bed aggradation processes prevailed, causing a total net deposition of 1016 m3.
This initial period of field research was chronologically followed by another type B PECP (September 2019 – July 2020), characterized by a lower entrainment capacity and a higher net bed excavation rate in both RCRs. This is corroborated by the data in
Table 9 and
Table 10, where the ANTD values vary between -0.16 and -0.19 m, and the net volume difference per 100 m2 (NVD*) is -19.3 m3 in UPR and -15.8 m3 in MDR. The incision experienced in this PECP was quickly counteracted in the following stage (July 2020 – February 2023), during which a class 2 event and a type D pattern caused the opposite effect: bed re-elevation due to a greater accumulation of bedload, which ended up exceeding the rate of transitory erosion. Nevertheless, the mean net deposition was much lower than that recorded in the first monitoring PECP, and the differences between the upper and middle reaches were also more noticeable (ANTD of 0.145 m in UPR and 0.026 m in MDR) (
Table 9 and
Table 10).
The other 12 PECPs were simulated, resulting in five PECPs with an ECP type B, two with ECP type C, four representing an ECP type D, and one of ECP type E. Retrospective analysis of these stages revealed a change in the trend of the erosion-deposition balance. From 1997 to 2010 a progressive lowering of the bed took place, interrupted only by a slight growth and stabilization along the middle reach during the events of November 19, 2003 and November 24, 2007, both of class 2. Between 2010 and 2013 there were three type D PECPs, each consisting of a single class 2 event, during which vertical bed accretion occurred in both stretches (NVD* of 15.7 m3 per 100 m2 in UPR and 7.3 m3 per 100 m2 in MDR). Finally, the 2013-2018 simulation period began with a PECP type B and ended with two class 4 minor events, which contributed to the generation of incision (NVD* = -26.7 m3 per 100 m2 in UPR and -15.9 m3 per 100 m2 in MDR). However, if we look at the frequency of PECPs with a positive or negative balance and the magnitude of morphological adjustment per stage, it is obvious that from the beginning of the simulated period there is a progressive change in trend, moving from linear erosion in the first nine PECPs (1997 to 2012) to the prevalence of deposition in the last six (2012 to 2022). The comparison of the averaged NDV* values for each interval of three consecutive PECPs confirms this in the two reference stretches, although in UPR the magnitude change is more pronounced. In this section, the average NDV* per 3 years changes from around -17 m3 per 100 m2 in the stage 1997-2007 to -6.3 m3 per 100 m2 in the 2007-2012 stage, and to 2.1 and 5.8 m3 per 100 m2 in the two most recent simulated PECPs triperiods (2012-2018 and 2018-2022, respectively).
Table 10.
Average and cumulative net difference in volume per 100 m2 and thickness difference for the changing surface area along the Upper and Middle RCRs.
Table 10.
Average and cumulative net difference in volume per 100 m2 and thickness difference for the changing surface area along the Upper and Middle RCRs.
|
UPR |
MDR |
PECP |
TAI |
TNVD |
NVD* |
CNVD |
ANTD |
CNTD |
TAI |
TNVD |
NVD* |
CNV |
ANTD |
CNTD |
July 2020 – Feb. 2023 |
2730 |
395.2 |
14.5 |
14.5 |
0.15 |
0.15 |
4532 |
118.7 |
2.6 |
2.6 |
0.03 |
0.03 |
Sep. 2019 - July 2020 |
2976 |
-575.3 |
-19.3 |
-4.8 |
-0.19 |
-0.04 |
4707 |
-742.2 |
-15.8 |
-13.1 |
-0.16 |
-0.13 |
Sep. 2018- Sep.2019 |
2763 |
613.2 |
22.2 |
17.4 |
0.22 |
0.18 |
4885 |
1013.1 |
20.7 |
7.6 |
0.21 |
0.08 |
Oct. 2013 – Sep 2018 |
2551 |
-680.5 |
-26.7 |
-9.3 |
-0.27 |
-0.09 |
4490 |
-716.1 |
-15.9 |
-8.3 |
-0.16 |
-0.08 |
Sep. 2012 – Oct. 2013 |
2543 |
443.8 |
17.5 |
8.2 |
0.17 |
0.08 |
4155 |
336.3 |
8.0 |
-0.3 |
0.08 |
0.00 |
Aug. 2012 – Sep. 2012 |
2543 |
395.3 |
15.5 |
23.7 |
0.16 |
0.24 |
4155 |
298.0 |
7.2 |
6.9 |
0.07 |
0.07 |
June 2010 – Aug. 2012 |
2543 |
356.8 |
14.0 |
37.7 |
0.14 |
0.38 |
4155 |
270.7 |
6.5 |
13.4 |
0.07 |
0.14 |
March 2009 – June 2010 |
2542 |
-755.9 |
-29.7 |
8.0 |
-0.30 |
0.08 |
4155 |
-894.3 |
-21.5 |
-8.1 |
-0.22 |
-0.08 |
Oct. 2007 – March 2009 |
2542 |
-78.5 |
-3.1 |
4.9 |
-0.03 |
0.05 |
4155 |
112.4 |
2.7 |
-5.4 |
0.03 |
-0.05 |
Nov. 2005 – Oct. 2007 |
2542 |
-802.9 |
-31.6 |
-26.7 |
-0.32 |
-0.27 |
4155 |
-983.7 |
-23.7 |
-29.1 |
-0.24 |
-0.29 |
Nov. 2003 – Nov. 2005 |
2542 |
-64.2 |
-2.5 |
-29.2 |
-0.03 |
-0.30 |
4155 |
85.0 |
2.1 |
-27.0 |
0.02 |
-0.27 |
March 2002 – Nov. 2003 |
2533 |
-496.9 |
-19.6 |
-48.8 |
-0.20 |
-0.50 |
4153 |
-771.8 |
-18.6 |
-45.6 |
-0.19 |
-0.46 |
Dec. 1998 – March 2002 |
2542 |
-847.4 |
-33.3 |
-82.1 |
-0.33 |
-0.83 |
4155 |
-987.5 |
-23.8 |
-69.4 |
-0.24 |
-0.70 |
Oct. 1997 – Dec. 1998 |
2533 |
-754.6 |
-29.9 |
-112.0 |
-0.30 |
-1.13 |
4153 |
-879.2 |
-21.2 |
-90.6 |
-0.21 |
-0.91 |
Oct. 1996 – Oct. 1997 |
2533 |
421.5 |
16.6 |
-95.4 |
0.17 |
-0.96 |
4153 |
319.2 |
7.7 |
-82.9 |
0.08 |
-0.83 |
The highest bed lowering rate was recorded in the PECPs from October 1997 to December 1998, December 1998 to March 2002, November 2005 to October 2007, and March 2009 to June 2010, with NVD* values of -29.7 to -33.3 m3 per 100 m2 and ANTD values of -0.30 to -0.33 m for UPR; and somewhat lower for MDR (NVD*: -21.2 to -23.8 m3 per 100 m2; ANTD: -0.21 to -0.24 m). In these stages, intense flash floods caused more incision than deposition. The comparison of the total volumes of surface lowering (VSL) with those of surface raising (VSR) during such sub-periods (
Figure 8) also corroborates this. By contrast, from 2012 the PCEPs included more moderate events, interspersed with bankfull and sub-bankfull peak discharges, which had high competence and effectiveness in transport. The bedload deposited during the recess stage of the hydrograph, in these cases, was generally greater than the sediment removal in the phase of transitory erosion. The most notable example corresponds to the PECP September 2018-September 2019, which was represented by a type A pattern and a single bankfull event. In this PECP, the highest vertical accretion rate of the entire analyzed period was recorded (NVD* = 22.2 m3 per 100 m2 and VSR = 613 m3 in UPR; NVD* = 20.7 m3 per 100 m2 and VSR = 1013 m3 in MDR) (
Figure 8).
This change in trend can be better appreciated by analyzing the variation in the cumulative net thickness differences (CNTD) and the cumulative unit volume of surface raising and lowering (CUVR and CUVL). The CNTD variation denotes a dominant bed downcutting process from October 1997 to June 2010, going from -1.13 m to 0.08 m in UPR and from -0.91 m to -0.08 m in MDR (
Figure 8). From June 2010 to February 2023, the bed experienced slight changes in elevation, maintaining a metastable dynamic equilibrium and a slight increase in net deposition along both stretches. The unit volumes of incision or aggradation do not take into account the bedform size or the affected area, but they do make it possible to detect the rate of morphological adjustment per surface unit and therefore to compare the sediment budgets in sections of different length. The most conclusive result, in this sense, is the verification that the rate of increase of CUVR was regularly uniform (0.8 m3 m-2 year-1), while the rate of incision followed an uneven evolution throughout the entire period, except in the 2010-2013 stage, in which it showed some similarity to the gradual bed elevation.
Figure 8.
The upper histograms represent the total volumes of lowered (VSL) and raised (VSR) surfaces for all PECPs; at the bottom, the left-hand figure shows the cumulative net thickness difference (CNTD), and the right-hand figure shows the profiles depicting the cumulative unit volume (CUV) of surface lowering (SL) and raising (SR) for the same sub-periods and channel reaches. AV = average values for each sub-period.
Figure 8.
The upper histograms represent the total volumes of lowered (VSL) and raised (VSR) surfaces for all PECPs; at the bottom, the left-hand figure shows the cumulative net thickness difference (CNTD), and the right-hand figure shows the profiles depicting the cumulative unit volume (CUV) of surface lowering (SL) and raising (SR) for the same sub-periods and channel reaches. AV = average values for each sub-period.
4.3. Retrospective Changes in Bed Elevation and Sediment Budgets Derived from TLS-Generated HRDEMs and the Resulting DoDs
The retrospective reconstruction of the morphological channel adjustments, performed at the event scale from TLS, provided more detailed information on topographic changes and bedload fluxes. Precise changes in the shape and extent of the bedforms can be seen in
Figure 9. Depending on the section and stage in question, such changes reflect different geomorphic dynamics. If we analyze the cumulative differences from 2020 to the past, in both PBSAs, during the three-year period (2017-2020) bed homogenization occurred, slightly prevailing deposition over erosion. As we go back to 2007 the topographic differences between pools and riffles in the upper PBSA increased, leaving bedform units more prominent. Effects accumulated around 2002 showed a quite braided pattern, as a result of the disaggregation of oblique and mid-channel bars and the formation of new active bars, associated with a high bedload. The incision of secondary channels in these older stages extended until 1997, but the longitudinal and median alluvial bars maintained their main features and only suffered small local changes.
Upstream, pool-riffle sequences appeared at all stages as this is either a supply- or transport-limited system with a gravel bed and a slope of around 0.02, which agrees with what was argued by [
51], [
52], and [
53]. Perhaps, the major morphological bed changes corresponded to narrow and non-sinuous channel stretches subject to geological control. Along these stretches such adjustments were promoted by high stream power values concentrated in short distances, as already suggested by Conesa-García et al. [
24] for entrenched channel cross-sections, along stretches with sudden changes in bed roughness.
The evolution of the MDR was always closely related to the geomorphic activity of the UPR and the intermediate section. Here again, previous incision features in secondary channels and the lateral edge of a mid-channel bar were blurred by deposits that filled them in the 2012-2016 sub-period. During this stage the bed remained practically stable as a result of an equilibrium in the balance between erosion and deposition. In recent years, it seems to be experiencing some vertical accretion and granular armoring that will most likely hinder the incision processes and favor lateral erosion.
Figure 9.
Cumulative HRDODs retrospectively simulated for the period 1997-2020 in the upper and middle RCRs, combining pre- and post-event ECs with TLS-derived DTMs over the monitoring stage (2018-2020).
Figure 9.
Cumulative HRDODs retrospectively simulated for the period 1997-2020 in the upper and middle RCRs, combining pre- and post-event ECs with TLS-derived DTMs over the monitoring stage (2018-2020).
Regarding the variations in bed elevation throughout the study period in UPR, aggradation was quite generalized along this stretch (
Table 11), except in the high alluvial bars that presented punctual erosion or are stabilized by large blocks and vegetation. Throughout this period (1996-2020) the pools were filled to a significant extent by sediments, but only locally did they reach the height of the riffles. The most noticeable filling took place in the 2007-2012 stage, with increases in bed elevation above 0.5 m in the deepest pools (
Figure 10). At the end of this stage, the bed reached its maximum topographic and morphosedimentary uniformity. In the last ten years, various events of different classes have occurred, causing uneven geomorphological impacts that have returned a certain irregularity to the bed. Several bars have been destroyed or displaced, and new pools have arisen, which denotes great geomorphic activity.
Table 11.
Average and cumulative changes in bed elevation relative to the current surface level (m), calculated by subtracting HRDoDs from TLS over simulated five-year periods.
Table 11.
Average and cumulative changes in bed elevation relative to the current surface level (m), calculated by subtracting HRDoDs from TLS over simulated five-year periods.
RCR |
Upper RCR |
Middle RCR |
Period |
ANTD |
CNTD |
St. Dev. |
ANTD |
CNTD |
St. Dev. |
2016-2020 |
0.06 |
0.06 |
0.04 |
0.04 |
0.04 |
0.06 |
2012-2016 |
-0.04 |
0.02 |
0.12 |
-0.01 |
0.03 |
0.09 |
2007-2012 |
-0.28 |
-0.26 |
0.33 |
-0.13 |
-0.10 |
0.21 |
2002-2007 |
-0.22 |
-0.48 |
0.58 |
-0.12 |
-0.22 |
0.37 |
1996-2002 |
-0.25 |
-0.73 |
0.82 |
-0.29 |
-0.51 |
0.63 |
Figure 10.
Serial profiles based on TLS data, showing the average thickness of the bed elevation differences along the channel center line over five or four-year sub-periods from 1996 to 2020. The area profiles represent the incision and aggradation sites for the entire period.
Figure 10.
Serial profiles based on TLS data, showing the average thickness of the bed elevation differences along the channel center line over five or four-year sub-periods from 1996 to 2020. The area profiles represent the incision and aggradation sites for the entire period.
The topographic changes were somewhat more moderate in MDR, although the events of each stage produced different specific effects. In the 2002-2007 stage, the erosion of the bar platforms and the accumulation in the pools and secondary channels were the most recurrent processes. Between 2007 and 2012 the erosion and mobilization of bedload material flattened and lowered the bed in the runs to elevation +0.2 m, while in the pools the surface rose +0.5 m compared to the previous stage. From 2012 to 2016 the bed remained relatively stable, with hardly any variation in its height, while in the most recent stage (2016-2020) vertical accretion at the tails of high bars, in pool sites, and at low bars occurred. Considering the entire period (1996-2020), it is worth noting the aggradation of the pre-existing bars in the growing development phase downstream (net deposition of 0.2 to 0.6 m) and the slight lowering of the supraplatform and head of the highest bars.
4.4. Spatio-Temporal Patterns in Stream Power Derived from Monitored and Simulated Peak Discharges
In
Table 12 the average values of the main energy indicators (ω, ωc, ∂ω/∂s, ε, and εc) are shown for each PECP and RCR. As already noted in section 3.5, these averages were estimated from datasets of 10 flow cross-sections in UPR and 12 along MDR, with a separation interval of 25 m. This type of data provides us with a simplified comparative vision of the interannual variation affecting the stream power in both sections, but, nevertheless, it masks its spatial variability. Detailed information on the spatial patterns of stream power in relation to channel changes in gravel-bed streams can be found in the works of Lea & Legleiter (2015) and Conesa-García et al. (2020, 2022). The latter authors found a wide range of ω and ωc values for our study stretches, which justified a large mosaic of spatial changes in bed elevations and sediment budgets attributed to events of different magnitude. Such variations were associated with substantial fluctuations in velocity and critical stress along the reference channel reaches, in accordance with the pattern suggested by [
51]. The spatial differences in roughness due to sudden changes in grain size and bedform also strongly influenced the calculated excess energy per unit bed area after specific events. A previous study on changes in stream power and morphological adjustments at the event-scale along these RCRs [
5] found larger differences between the standard deviations and the means of ω, ∂ω/∂s, and ε for Class 2 and 3 events than for large flash floods. Also, in MDR residual stream power values were obtained further from the mean than in UPR. All of this leads us to consider with more reserve the averaged values for PECPs that do not include overtopping discharges but do contain bankfull and sub-bankfull type events in their sequence (Pattern B). Taking this into account,
Table 12 indicates that the PECPs with the highest mean values of ω occurred between October 1997 and June 2010, highlighting the stages December 1998 – March 2002 and November 2005 – October 2007. In both, the average stream power exceeded 265 and 240 W m-2 in UPR and MDR, respectively.
Table 12.
Average values of stream power variables estimated for each PECP sub-period, for the Upper and Middle RCRs.
Table 12.
Average values of stream power variables estimated for each PECP sub-period, for the Upper and Middle RCRs.
PECP |
ω |
δω/δs |
ε |
εc |
UPR |
MDR |
UPR |
MDR |
UPR |
MDR |
UPR |
MDR |
JUL20-FEB23 |
204.9 |
122.2 |
1.01 |
-0.76 |
89.40 |
73.45 |
1.34 |
0.73 |
SEP19-JUL20 |
156.6 |
256.1 |
0.38 |
-0.79 |
72.35 |
59.42 |
0.84 |
1.19 |
SEP18-SEP19 |
214.1 |
154.2 |
0.82 |
-1.06 |
108.15 |
50.22 |
1.25 |
0.79 |
OCT13-SEP18 |
212.1 |
129.8 |
1.22 |
-1.35 |
98.92 |
39.15 |
1.44 |
0.69 |
SEP12-OCT13 |
117.1 |
270.4 |
0.64 |
-1.01 |
47.91 |
85.33 |
1.05 |
1.37 |
AUG12-SEP12 |
210.2 |
104.1 |
1.08 |
-1.42 |
92.50 |
31.02 |
1.57 |
0.62 |
JUN10-AUG12 |
184.4 |
113.1 |
0.60 |
-0.51 |
91.01 |
33.63 |
1.05 |
0.66 |
MAR09-JUN10 |
218.4 |
148.3 |
0.78 |
-1.62 |
104.01 |
49.31 |
1.20 |
0.76 |
OCT07-MAR09 |
230.8 |
183.8 |
1.24 |
-1.21 |
133.17 |
53.96 |
1.59 |
0.86 |
NOV05-OCT07 |
265.5 |
240.2 |
1.06 |
-1.21 |
126.97 |
82.78 |
1.47 |
1.09 |
NOV03-NOV05 |
153.6 |
112.9 |
0.38 |
-0.37 |
73.31 |
77.56 |
0.85 |
1.03 |
MAR02 -NOV03 |
143.6 |
102.8 |
0.34 |
-0.20 |
66.12 |
69.44 |
0.76 |
0.92 |
DEC98-MAR02 |
270.7 |
240.4 |
1.34 |
-1.52 |
157.30 |
120.86 |
1.82 |
1.12 |
OCT97-DEC98 |
261.7 |
234.2 |
0.82 |
-0.95 |
117.06 |
70.75 |
1.15 |
0.95 |
OCT96-OCT97 |
144.5 |
97.1 |
0.26 |
-0.15 |
64.46 |
48.38 |
0.75 |
0.74 |
Regarding the mean stream power gradient (∂ω/∂s), the contrast between the positive values of UPR and the negative values found in MDR is noteworthy. This is consistent with the assumptions described by Dade and Friend [
54] and Gartner et al. [
55] when analyzing the influence of stream power gradients on downstream sediment flux during floods. Downstream decreases in the average stream power gradient (negative ∂ω/∂s) were normally associated with depositional responses to peak discharges in the recession branch of their hydrograph. It is also of note that, depending on the pattern of event classes, the contrast of the mean values of ∂ω/∂s between the upper and middle reaches was greater or lesser. The greatest negative gradients averaged in MDR (∂ω/∂s < -1.20 W m-2 m-1) did not always correspond to the greatest positive ones obtained in UPR (∂ω/∂s > 1.20 W m-2 m-1). They had a close relationship with the latter in the stages with ECPs of type B and E, devoid of Class 1 events (e.g. October 2007 – March 2009, December 1998 – March 2002), but, on the other hand, they also coincided with lower ∂ω /∂s values in UPR (∂ω/∂s of 0.78 to 1.1 W m-2 m-1) during PECPs of type A, C, and D (e.g. September 2018 – September 2019, March 2009 – June 2010, August 2012 – September 2012).
The average values of ε and εc per PCEP were, as expected, higher in UPR than in MDR, except in three cases (March 2002 – November 2003, November 2003 – November 2005, and September 2012 - October 2013). Although bed roughness and critical stream power had a great influence on both variables in most of the retrospectively simulated periods, both ε and εc followed the ω variation patterns (
Figure 11), reaching their maximums in UPR during the PECPs December 1998 – March 2002, November 2005 – October 2007, and October 2007 – March 2009 (ε of 125 to 160 W m-2 and εc of 1.45 to 1.80 MJ).
Figure 11.
ω-ε and ω-εc relationships (left-hand graphs) and ∂ω/∂s-ε and ∂ω/∂s-εc relationships (right-hand graphs), at a significance level (p-value) < 0.05. AVP = average values per PCEP sub-period.
Figure 11.
ω-ε and ω-εc relationships (left-hand graphs) and ∂ω/∂s-ε and ∂ω/∂s-εc relationships (right-hand graphs), at a significance level (p-value) < 0.05. AVP = average values per PCEP sub-period.
4.5. Relationships between Morphological Bed Adjustments and Temporal Changes in Potential Sediment Inputs and Stream Power
The average (ASYP) and accumulated (CSYP) sediment yield per sub-period were plotted against ANTD, UVSL, and UVSR (
Figure 12) to examine the hypothesis that these variables could be correlated. All these scatter plots display point clouds, which represent such types of relationships for the PECPs, including events monitored over the period 2018-2022 and those predicted retrospectively by GeoWEPP up to 1996. In general, the scatter diagrams show a large dispersion of residual values and there appears to be no correlation between morphological bed adjustment classes and sediment production at the basin level. Only in MDR is there a slight fit by polynomial regression of order 2 (r2 ≈ 0.65), that improves significantly at order 4 (r2 ≈ 0.87), when relating the average sediment yield to ANTD and UVSL. The other relationships present scattered and skewed distributions, with groupings of points that, relatively, are somewhat more organized for the cases in which the variable ASYP is considered. In fact, it is clear that ANTD has a chaotic relationship with CSYP, but shows a positive correlation with ASYP. The total lack of relationship between CSYP and ANTD is explained if it is taken into account that the net variations in bed elevation are due more to the effects of specific chronologically consecutive events of different classes than to the cumulative production of sediments in a given period.
Regarding the relationships between the unit volumes of net erosion or deposition and the average or accumulated values of sediment delivery, no clear trends of unique statistical populations are observed, but rather dual behaviors in both sections. Only in the middle reach can a decrease in UVSL as ASYP increases be verified. In the other cases the same pattern seems to be repeated: two groups per RCR with unequal behavior in relation to different contributions of sediments from their respective catchment areas. For both groups of PECPs, bed lowering occurred in UPR with less sediment generated at the basin level than in MDR. The first group identified in UPR includes events that originated a moderate mean sediment production (1000 to 3100 tn) and an average unit scouring volume of 0.22 to 0.35 m3 m-2. This same group of PECPs exhibited somewhat lower bed downcutting rates in MDR (0.23 to 0.28 m3 m-2), associated with slightly higher ASYP values (2000 to 4100 t). The second set of points reflected UVSL growth with ASYP, which generally showed a positive relationship through a gently sloping trend line, representing a low range of bed incision, accompanied by lower average sediment supply in UPR than in MDR.
Figure 12.
Average and cumulative sediment yield per PECP sub-period (ASYP and CSYP, respectively), estimated at the entrance of the upper (UPR) and middle (MDR) channel reaches, versus ANTD, UVSL, and UVSR.
Figure 12.
Average and cumulative sediment yield per PECP sub-period (ASYP and CSYP, respectively), estimated at the entrance of the upper (UPR) and middle (MDR) channel reaches, versus ANTD, UVSL, and UVSR.
A different pattern is observed when relating the average rates of vertical bed accretion to those of sediment production in the watersheds of each channel reach. In UPR, the points are concentrated in a small cloud in the lower left part of the scatter plot, which translates into a succession of PECPs with low average net deposition (UVSR of 0.1 to 0.23 m), associated with little productive activity in the sediment source areas (ASYP < 3500 t). By contrast, MDR experienced greater variability in the rate of bed rise (UVSR of 0.13 to 0.37 m) and greater availability of sediments for transport and deposition (2100 < ASY < 9800 t). As a result, this RCR frequently acted as a sink for material eroded from the closest upstream reach. This is consistent with the results obtained at a large scale by Wilcock and Crowe [
56] and Török et al. [
57], using flume experiments with mixed-size bed sediment, in which eroded particles tended to be deposited immediately after the erosion zone. The same relationship patterns described for ASYP are repeated for CSYP, although here the groupings are linear and show a certain positive trend in relation to the lowering and growth of the bed.
The point clouds in
Figure 13 show the ANTD values (average, standard deviation, and unit error) in relation to the mean stream power gradient (∂ω/∂s), and to the excess energy per unit bed area (ε) expended above ωc and accumulated (εc), averaged for each of the periods defined by PECPs. In the first column, positive and negative values of ANTD depict surface lowering and raising, respectively. Scatter plots relating ANTD to the mean ∂ω/∂s ratio represent two very different trends in each stretch: a positive relationship in UPR and a negative one in MDR, with a separation threshold at ∂ω/∂s = 0. In both RCRs, the closer the ∂ω/∂s values are to this threshold, the greater the erosion that is observed (ANTD < -0.3 m in UPR and ≈ -0.2 m in MDR); in contrast, extreme values, both negative and positive, correspond to the highest rates of net deposition (ANTD > 0.2 m). The field surveys carried out during the monitoring period already showed morphological adjustments in both directions, deposition and bed scouring, depending on the magnitude of each event and the availability of sediments for transport. After the largest flood in this period (19 April, 2019), the upper reach experienced greater net erosion than the middle reach, most likely especially in the sections where ω increased downstream above the ωc, coinciding with positive ∂ω/∂s values. This specific case was described by Conesa-García et al. [
5] when explaining the effect of the energy gradient on the net flux variations immediately downstream.
Figure 13.
Values of stream power variables (∂ω/∂s, ε, and εc) versus parameters related to the net thickness difference (ANTD, SD-NTD, and Er-NTD). AVP = average values per PCEP sub-period; AOV and ASV are, respectively, the observed and simulated mean values for each sub-period. Er-NTD expressed per unit (p.u.).
Figure 13.
Values of stream power variables (∂ω/∂s, ε, and εc) versus parameters related to the net thickness difference (ANTD, SD-NTD, and Er-NTD). AVP = average values per PCEP sub-period; AOV and ASV are, respectively, the observed and simulated mean values for each sub-period. Er-NTD expressed per unit (p.u.).
Deposition was slightly greater in the middle stretch than in the upper reach, but the point distributions relating ANTD to ∂ω/∂s, ε, and ε
c differed considerably between the two channel reaches. Among all the relationship patterns of these variables, the most clearly defined are those relating ANTD to ∂ω/∂s, in both RCRs. In this case, UPR and MDR show different trends supported by linear regressions with good fits (r
2 of 0.64 and 0.72, respectively, with a significance level, p-value, < 0.05). For UPR there is a positive correlation between ANTD and ε
c (r
2 = 0.65) (
Table 13). The remaining relationships have poorly fitted distributions, although in certain cases they show specific behavior patterns depending on the type of section. Both UPR and MDR present two groups of residual values, one positive and the other negative, with a slightly progressive increasing trend. However, some differences can be noted. For example, in the middle section, sedimentary vertical accretion began to occur at a much lower excess energy than in the upper stretch: in UDR, a mean ε of 90 Wm
-2 and an ε
c of 1.2 MJ were needed for net deposition to occur in a PECP, while in MDR 30 Wm
-2 and 0.6 MJ respectively, were enough. Assuming that these relationships follow average behavior patterns for the simulation periods, heterogeneous distributions relating the energy variables to the dispersion (SD-NTD) and unit error (Er-NTD) values of ANTD are found here, which could explain the great variability of the underlying geomorphological processes within each period.
Table 13.
Significant relationships between bed morphological change parameters and stream power variables for the upper and middle RCRs (p-value < 0.05 for the confidence interval of 95%).
Table 13.
Significant relationships between bed morphological change parameters and stream power variables for the upper and middle RCRs (p-value < 0.05 for the confidence interval of 95%).
RCR |
Relationship |
Regression equation |
Function |
r2
|
UPR |
UVSL vs ∂ω/∂s |
UVSL = 0.093 ∂ω/∂s2 – 0.363 ∂ω/∂s + 0.419 |
Polynomial |
0.73 |
ANTD vs ∂ω/∂s |
ANTD = 0.485 ∂ω/∂s – 0.450 |
Linear |
0.64 |
ANTD vs εc
|
ANTD = 0.524 εc – 0.698 |
Linear |
0.65 |
UVSC vs ∂ω/∂s |
UVSC = 0.192 ∂ω/∂s2 - 0.533 ∂ω/∂s + 0.638 |
Polynomial |
0.75 |
UVSL vs εc |
UVSL = 0.178 εc2 – 0.670 εc + 0.732 |
Polynomial |
0.75 |
UVSC vs εc
|
UVSC = 0.319 εc2 – 1.027 εc + 1.10 |
Polynomial |
0.78 |
MDR |
ANTD vs ∂ω/∂s |
ANTD = -0.269 ∂ω/∂s – 0.309 |
Linear |
0.72 |
ANTD vs ASY |
ANTD = -5E-09 ASY2 + 0.0001 ASY – 0.48 |
Polynomial |
0.67 |
UVSL vs ∂ω/∂s |
UVSL = -0.076 ∂ω/∂s2 + 0.003 ∂ω/∂s + 0.282 |
Polynomial |
0.81 |
UVSR vs εc
|
UVSR = 0.608 εc2 – 0.905 εc + 0.479 |
Polynomial |
0.89 |
UVSL vs ASY |
UVSL = 121.01 ASY-0.78
|
Power |
0.64 |
The variation patterns of the unit rates of lowering (UVSL) and raising (UVSR) and total bed change (UVSC) in relation to the stream power indicators are characterized by a different level of dispersion and grouping, depending on the channel reach (
Figure 14). Here again, the best fits are found for the relationships between the bed modification unit volumes and the mean stream power gradient and the variables derived from the energy balance (
Table 13). These are polynomial regression fits, indicating a nonlinear relationship between the ∂ω/∂s values and the corresponding conditional means of UVSC (in UPR) and UVSL (in MDR), or between ε
c and the average UVSC (in UPR) and UVSR (in MDR).
Figure 14.
Values of stream power variables (∂ω/∂s, ε, and εc) versus unit volume of surface lowering (UVSL), raising (UVSR), and changing (UVSC). The values are averaged per PCEP sub-period (AVP).
Figure 14.
Values of stream power variables (∂ω/∂s, ε, and εc) versus unit volume of surface lowering (UVSL), raising (UVSR), and changing (UVSC). The values are averaged per PCEP sub-period (AVP).
Along the middle RCR, the UVSR values corresponded to negative ∂ω/∂s values, the highest being those closest to zero, a threshold around which the greatest deposition also occurred in UCR, decreasing as the ∂ω/∂s values became more positive. This bears some resemblance to the trend described by Lea and Legleiter [
58] for perennial gravel streams and by Conesa-García et al. [
5] at the event scale for the RCRs studied here. The worst fits, on the other hand, are those relating the unit morphological bed adjustments to the mean stream power, especially in the case of global changes that seem not to be justified solely by the flow competence, but also require consideration of the critical stress and bed roughness.
An interesting aspect can be appreciated in the scatter plots in
Figure 14 showing the relationships between the mean unit volume of surface lowering or raising and the average and accumulated values of excess energy for the PECP periods. The point clouds relating UVSL and UVSR to ε and εc are generally quite dispersed but despite this we can observe some differences. Specifically, two slightly different patterns are distinguished when comparing the trends of these relationships and value groups according to the RCR. In UPR, a progressive increase in εc seems to coincide with a decrease in the unit bed scour rate, while in MDR it implies an increase in UVSR. This has a certain consistency since the floodwaters in the upper reach are usually clearer and tend to cause bed incision.
On the other hand, downstream, in the middle section, the flow energy, which increases with the magnitude of the event, is used to transport a large bedload that ends up deposited when the velocity decreases. If, on the contrary, we consider the average ε values for each of the simulated stages in the last 27 years, dispersion is the most common feature of the distributions in both RCRs, perhaps as a result of discontinuous changes in bedload and channel morphology produced during this period. The scattered and skewed distributions of this variable could be related to nonhomogeneous bedforms or changes in granular texture. Conesa-García et al. [
24] associated this lack of relationship in another complex gravel-bed dry channel (upper Mula stream, in southeastern Spain) with the presence of blocks from the bank breaks, pools-riffle sequences, and local transitions from alluvium to substrate outcrop and vice versa. Zapico et al. [
59] analyzed this type of relationship in a steep, sand-gravel ephemeral channel and found, by contrast, a clear relationship between changes in bedload flux and morphological channel adjustments.
As a result of analyzing these relationship patterns, it can be stated that rather than linear equations they correspond to polynomial regressions derived from different geomorphological dynamics and processes. These are manifested on a temporal and spatial scale through a wide range of feedback effects involving external environmental factors and are intrinsic to the channel itself. In many of the patterns described, there does not seem to be a clear relationship between sediment production at the basin level and morphological channel adjustments. However, an inverse relationship has been found between CSY in each PECP and εc in the same period (
Figure 15).
The higher potential inputs of basin sediments recorded in MDR do not correspond to greater vertical bed accretion in all PECPs or to an increase in εc along this stretch. This is most likely due to a marked decrease in the ω/ωc ratio in the middle section with respect to the upper one, linked to the decrease in downstream slope and increase in bed roughness. However, the largest flood events always coincided with the highest peak flows (30 > Qp < 80 m3 s-1, class 1, overbank), the greatest accumulated amounts of sediment yield at the basin level (12 < CSY < 27 t), the highest values of εc inside the channel (1.8 < εc < 3.7 MJ), and the most notable bedform changes for both channel reaches (20 < NVD* ≈ 22 m3 per 100 m2 on 19 March, 2019; and NVD* of -21 to -33 m3 per 100 m2 after the three most important flash floods during the 2000-2009 period). Moderate and minor events, class 3 and 4, respectively, generated insignificant bed modifications related to low accumulated sediment production and εc values below 0.6 MJ. The upper stretch also showed peculiar behavior during the bankfull stages, when SY did not exceed 3000 t; but, εc was much higher than in MDR, reaching values above 1.3 MJ, and the net volume differences in bedforms were especially significant (NVD* ≈ -20 m3 per 100 m2 in 2003-2003 and 14 and 18 m3 per 100 m2 in 2012-13 and 2022-23, respectively).
Figure 15.
Comparison of the sediment yield and excess accumulated energy estimated for each of the events that occurred in the period 1996-2022 along UPR and MDR.
Figure 15.
Comparison of the sediment yield and excess accumulated energy estimated for each of the events that occurred in the period 1996-2022 along UPR and MDR.