Introduction

Rapid channel incision reflecting a loss of geomorphic dynamic equilibrium is a common phenomenon in contemporary fluvial systems affected by intense human disturbances (Florsheim et al. 2013). Incision leads to increased flow capacity of the channel and results from imbalance between sediment load and the river’s ability to transport that load, and from the occurrence of channel boundary conditions (such as bank reinforcements) that prevent re-establishment of a state of equilibrium via river metamorphosis (Wyżga et al. 2016c). River incision may thus result from a single causal factor, such as in-channel gravel mining (Peiry 1987; Collins and Dunne 1989), disruption of the continuity of bedload transport by a dam reservoir (Williams and Wolman 1984), channel regulation (Brookes 1987; Simon 1989) or land-use change in the catchment (Liébault and Piégay 2001; Keestra et al. 2005), or from a combination of different factors (e.g. Landon et al. 1998; Rinaldi 2003; Surian and Rinaldi 2003; Škarpich et al. 2013).

During the twentieth century, channel incision was pronounced in Polish Carpathian rivers (Klimek 1983; Wyżga 2008; Wyżga et al. 2016c). Large-scale mining of gravel from channels (Rinaldi et al. 2005; Wyżga et al. 2010), an increase in river transport capacity resulting from channelization works (Wyżga 2001a; Korpak 2007; Zawiejska and Wyżga 2010), and a decrease in sediment supply from the catchments subjected to land-use changes (Lach and Wyżga 2002) were indicated as the causes of rapid incision of these rivers.

Channel incision brings about a number of detrimental hydraulic, hydrological and ecological effects (Bravard et al. 1997; Wyżga 2008). A lowering of water stages may result in the exposure of bank water intakes above medium and low flows. Concentration of flood flows in the narrow channels of high flow capacity results in high values of flow velocity, unit stream power and bed shear stress (Radecki-Pawlik et al. 2016), which facilitate undermining of bank protection structures and bridge foundations. In-channel concentration of flood flows reduces floodwater retention in floodplain areas (Froehlich 1980; Wyżga 2008) and overbank sedimentation on valley floors (Wyżga 2001b, 2008). The reduction in floodwater retention on floodplains decreases attenuation of flood waves and leads to increased flood hazard in downstream river reaches (Wyżga 1997; Wyżga et al. 2016a). Disconnection of channels from their floodplains (Wohl 2004) disrupts lateral connectivity of river ecosystems, whereas transformation of the alluvial channel beds into bedrock ones—frequent in incised Polish Carpathian rivers—eliminates vertical connectivity of the ecosystems (Hajdukiewicz et al. 2019). Such changes deteriorate the ecological status of incised rivers (Bravard et al. 1997).

A large number of detrimental effects of channel incision call for action to arrest or reverse the negative changes of river channels. Implementation of effective remedial measures to stop and reverse the tendency of rivers to incise is feasible; however, it requires earlier recognition of the scale of channel incision along river courses. To date, the amount of channel incision of Polish Carpathian rivers has been determined only for water-gauge cross sections, where it was reconstructed on the basis of long-term lowering of minimum annual stages of the rivers (Wyżga 2008; Wyżga et al. 2016c; Hajdukiewicz et al. 2019). However, this degree of recognition of the phenomenon is inadequate for remedial actions because neighbouring river sections may be typified by different amounts of accomplished incision or even by the occurrence of considerable incision in one section and its absence in the other (Zawiejska and Wyżga 2010). This is especially true with respect to mountain river reaches, where the exposure of resistant bedrock layers on the channel bed or the construction of weirs may result in radical differences in the amount of channel incision between adjacent river sections.

A survey method allowing for recognition of the amount of channel incision along a river is photogrammetry. It was successfully used to assess changes in the relief of sand dunes in the Słowiński National Park, Poland (Michalowska and Glowienka 2008). This study used five digital elevation models (DEMs) extracted for the dates when photogrammetric flights were performed, as well as differential DEMs. The detection of changes in dune relief on the basis of photogrammetric survey was relatively easy because of substantial changes in vertical position of dunes between the dates and the lack of vegetation. In another study, a time series of DEMs was also analysed to examine the response of long profiles of watercourses flowing to the Dead Sea, Israel, to rapid and continuous base-level lowering (Bowman et al. 2010). Application of the method documented high rates of incision of the channels into exposed cohesive sediments of the Dead Sea during the last three decades. DEMs created from sequential aerial photos were also used to reconstruct gully evolution on bare, deforested hillslopes in New Zealand (Derose et al. 1998; Betts et al. 2003) and on moraines recently exposed after glacial retreat in British Columbia, Canada (Schiefer and Gilbert 2007), and to assess the amount of sediment redistributed as a result of the gully development.

Environmental conditions typifying Polish Carpathian rivers differ from those in the areas considered in the above-mentioned studies. The presence of vegetation and relatively small vertical changes between the dates of consecutive aerial photos necessitate a higher accuracy of photogrammetric survey. To verify the usefulness of photogrammetric surveys for identification of the amount of changes in vertical position of a river channel over the last few decades, we analysed such changes in a 3-km-long reach of the Czarny Dunajec, for which a few sets of archival aerial photos from the period 1964–1994 together with recent orthophotos and DEMs were available. In this reach, the river has experienced considerable but spatially variable changes in the vertical position of its channel since the mid-twentieth century. In selected cross sections from this reach, vertical position of paleochannels from different periods of the twentieth century was previously determined by means of ground surveys (Zawiejska and Wyżga 2010; Wyżga et al. 2016b) and this allows for approximate verification of the photogrammetric reconstruction of changes in vertical river position.

Study area

The Czarny Dunajec (Fig. 1) is the upper part of the Dunajec, one of the largest rivers in the Polish Carpathians. It rises at ~ 1500 m a.s.l. in the Tatras, a high-mountain massif in the Western Carpathians. The high-mountain part of the catchment (with the highest peak at an altitude of 2176 m) determines the hydrologic regime of the river characterized by low flows between November and March, a prolonged but relatively low snowmelt freshet in late April and May, and highly variable flows during the warm part of the year, with a possible occurrence of a few days-long, relatively high floods between May and August. In the Tatra Mountains foreland, the river flows for 38 km, first within the Gubałówka Hills and then through the Orawa–Nowy Targ Basin. In this reach, the river flows within a non-cohesive alluvial plain formed from quartzitic, granitic and sandstone material, with an admixture of calcareous particles along first 10 km downstream of the Tatras.

Fig. 1
figure 1

Drainage network of the Czarny Dunajec catchment and detailed setting of the study reach (a). Location of the Czarny Dunajec River in relation to physiogeographic regions of southern Poland (b) and in the country (c). 1—high mountains; 2—mountains of intermediate and low height; 3—foothills; 4—intramontane and submontane basins; 5—boundary of the Czarny Dunajec catchment; 6—boundaries of physiogeographic units; 7—flow-gauging stations; 8—study reach

Until the mid-twentieth century, the river in its foreland reach flowed in a wide, shallow, multi-thread channel (Zawiejska and Wyżga 2010). In the second half of the century, the river experienced a number of disturbances (Zawiejska and Krzemień 2004; Zawiejska and Wyżga 2010; Wyżga et al. 2012). In the 1950s–1960s, large-scale gravel mining was conducted at a few locations, with the activity focused on the extraction of larger particles. The intensity of material removal from the river can be evidenced by the construction of rail tracks at three locations to transport the material out of the channel by train (Dudziak 1965). In the following decades, this large-scale mining was replaced by illegal, common extraction of cobbles and boulders from the channel at various locations (Krzemień 2003; Wyżga et al. 2010). Even though this last activity removed a relatively small amount of material from the river, it also exposed to flow finer grains from the subarmour layer, which could be then relatively easily flushed out to downstream river sections (Zawiejska et al. 2015). Channelization works carried out along most of the river course resulted in the replacement of the wide, multi-thread channel by a narrow, single-thread one (Krzemień 1981; Zawiejska and Wyżga 2010; Hajdukiewicz and Wyżga 2019). This increased transport capacity of the river and stopped sediment delivery from eroded channel banks. Finally, cessation of pasture and forest grazing in the Tatras following the establishment of a national park in the high-mountain part of the catchment in 1955 and a decrease in the pastoral and agricultural use of hillslopes in its foreland part reduced catchment sediment supply to the river (Wyżga et al. 2012). A combined response of the river to these impacts was rapid channel incision.

A 3-km-long river reach between the villages of Chochołów and Koniówka (Figs. 1a, 2) is typified by especially deep channel incision. This reflects intense gravel mining carried out close to the downstream end of the reach in the 1950s–1960s and lithological characteristics of the bedrock. Only in the uppermost part of the reach located within the Gubałówka Hills, the river is underlain by relatively resistant sandstone layers of Paleogene age. Along most of the reach, located within the intramontane Orawa–Nowy Targ Basin, the river is underlain by unconsolidated clays of Miocene–Pliocene age. After the dissection of gravelly alluvium, the river continued to incise into the clays, while maintaining a thin layer of gravels on the channel bed (Zawiejska and Wyżga 2010). As rapid channel incision became apparent shortly after the onset of large-scale gravel mining, in the early 1960s an 1.8-m-high weir was constructed at the downstream end of the reach to prevent undermining of the piers of the bridge at Koniówka. In the subsequent decades, vertical position of the river channel in a 7-km-long downstream reach was stabilized by a number of concrete weirs. In turn, immediately upstream of the Chochołów–Koniówka reach, the river incised by about 1 m, after which bed degradation ceased with the exposure of resistant sandstone layers on the channel bed.

Fig. 2
figure 2

Study reach of the Czarny Dunajec shown on orthophotos from 1964 and 2009

During the second half of the twentieth century, forest cover on the floodplain in the Chochołów–Koniówka reach gradually increased from 40% in 1964 to 87% in 1994, and since then, it remained practically unchanged (Hajdukiewicz and Wyżga 2019). Hydrological data recorded since 1968 at the Koniówka water-gauge station (Fig. 1a) characterize conditions in the reach. Here, the average for the lowest annual discharges is 0.91 m3 s−1, mean annual discharge equals 4.4 m3 s−1, and the average for the highest annual discharges is 56 m3 s−1. Notably, despite considerable inter-decadal variability in flood magnitudes on the river (Ruiz-Villanueva et al. 2016), no significant trend during the observation period was found for maximum annual discharges (R2 = 0.001; p = 0.85), mean annual discharges (R2 = 0.002; p = 0.76) and minimum annual discharges (R2 = 0.01; p = 0.46) recorded at the station.

Materials and methods

Photogrammetric materials

Four sets of archival aerial photos were available for the study reach—one for each decade between the 1960s and the 1990s (Table 1). Three of them were panchromatic and one (from 1994) in colour. Due to the scale of the photos, the whole study reach was covered by DEMs extracted from a single stereopair of the photos from 1964, 1983 and 1994. The photos from the neighbouring rows were used for improvement of external orientation of aerial photos and generating control DEMs. While the altitudes of the flights varied, the parameters of the cameras were comparable. Standards for automatic extraction of a DEM from aerial photos are established, depending on terrain relief, as 0.25–10% of a flight altitude above the ground level (GUGiK 2001). Thus, for the moderately flat valley floor of the Czarny Dunajec River, the expected accuracy of DEM extraction from the used photos varies between 0.25 and 0.45 m. The terrain pixel dimensions vary between 15 and 42 cm depending on the photo scale and on scanning resolution.

Table 1 Parameters of archival aerial photos used in the study

The reference materials for the assessment of changes in the vertical position of the river channel consisted of orthophotos from 2009 and 2012 and DEMs: stereophotogrammetrically derived one from 2009 and ALS-derived from 2012. These materials satisfy high-accuracy requirements for topographic data used in flood hazard assessment in Poland (Kurczyński and Bakuła 2013). High accuracy of ALS-derived DEM for forested areas was also proven in field conditions (Hajdukiewicz and Romanyshyn 2017).

Photogrammetric processing

At the stage of aerial photo processing, a comparison of aerotriangulation conditions for sequential aerial photos was performed for the study reach and two additional test areas: one located ca. 10 km downstream of the studied reach of the Czarny Dunajec and one situated in the Białka River valley 30 km to the east (Hajdukiewicz 2015). The aim of this analysis was to determine different kinds of distortion (radial, tangential, and irregular, caused by storage conditions) of the aerial photos. For photogrammetric projects established for each of the areas and given dates when aerial photos when taken, a number of control points were surveyed using different methods. Spot heights measured in the field with the GNSS RTN method were used to assess the accuracy of each DEM generated from a set of archival aerial photos. Objects stable for decades, such as intersections of property boundaries, corners of buildings or bases of low voltage network poles were selected as control points. The final selection of the control points in the study area was based on the experiment made in the test areas. Objects impossible to identify in the field were used as check points for the DEMs in the second stage, after identified objects were used as ground control points (GCP) for performing external orientation of the aerial photos. Moreover, a number of ground control points measured on a vector map and the orthophoto of the area were used. In such cases, horizontal coordinates of the GCPs were acquired from these source data, but vertical coordinates were read on the detailed ALS-acquired DEM from 2012. This combined method allows for sufficient accuracy of a GCPs survey (Michalowska and Glowienka 2008). It was assumed that the accuracy of extracted DEMs strongly depends on the accuracy of aerotriangulation, expressed by RMS error after bundle adjustment process. Therefore, a number of iterations of bundle adjustment were performed with PCI orthoengine software, using different sets of control points for each set of aerial photos. Finally, RMS errors did not exceed 0.25 m for any date of analysed aerial photos that meets requirements of photogrammetric surveys (GUGiK 1999). Generally, a strong dependence of the observed error on the age of aerial photos was found, caused by such factors as the radiometric quality of the photos, geometric distortion of the photo material (resulting from improper storage conditions) and uncertainty of location of ground control points.

Extracting and editing of DEMs

DEMs were extracted using PCI orthoengine software, with multiple iterative algorithm of matching on descending image overview levels. Accuracy of the DEMs was tested using the Point Sampling tool in Quantum GIS software. The mean error of ground observations on the DEM generated for test areas, referenced to ALS or GPS measurements, typically ranged between 0.3 and 0.4 m, which can be considered sufficient (GUGiK 2001). However, local maximum errors, probably caused by insufficient radiometric and geometric quality of the photos, exceeded 1 m. Spatial distribution of the errors was analysed using a set of 96 check points with bilinear interpolation of the plane of average error value. If a significant trend of the spatial distribution of errors was found (Fig. 3a), the next iteration of external orientation of the aerial photos was performed with a modified set of GCPs to obtain a new extracted DEM without a significant systematic error (Fig. 3b). After this procedure, local errors were considerably reduced, while remaining substantially larger than the overall trend. In this case, differences in the trend were comparable to the accuracy of the reference ALS survey.

Fig. 3
figure 3

Spatial distribution of errors on the DEM from 1964 extracted without known parameters of camera distortion (a) and after calculation of distortion parameters, with a modified set of GCP (b). The errors are shown on the background of orthoimage of a single aerial photo. 1—values of vertical error of check points (m); 2—isolines of vertical error; 3—study cross sections

Frequency of local errors depended on the date of photogrammetric material. They were rare on the DEM generated from the photos taken in 1977, more frequent on the DEMs extracted from the photos taken in 1964 and 1983, and the most numerous in the case of the DEM generated from the photos from 1994. In this last DEM, the errors coincided with grassland areas. DEM errors occurred also in forested areas because of improper measurement during the automatic extraction process where the surface of trees or shrubs must have been confused with the terrain surface (Fig. 4a). Such errors were corrected using the PCI Focus tool. The erroneously classified areas were masked using on-the-screen vector polygon digitization (Fig. 4b), and then, one of the following methods of altitude correction was applied for each polygon: linear interpolation from the edges in the case of riparian forest with variable tree heights or subtraction of a constant value in the case of trees or shrubs of the same height. The latter method was applied only to small areas of planted trees in the neighbourhood of riparian forest (ca 5% of edited area). Figure 4c presents the DEM after editing. The proportion of edited area equalled less than 1% for the DEM extracted from the photos taken in 1964, 2% for the DEM from 1983 and 6% for that from 1994. The increase reflected progressive expansion of riparian forest along the river (Hajdukiewicz and Wyżga 2019) and growing density of trees in the forest over the study period—in the case of sparse forest, DEM was usually extracted properly.

Fig. 4
figure 4

Fragment of the orthophoto of the Czarny Dunajec valley floor from 1994 (a) and DEM for this area after extraction of a stereopair of archival aerial photos (b) and after edition (c). Black lines delineate areas where the elevation of terrain surface was distorted by vegetation

River measurements on orthophotos and DEMs

A river centreline in the study reach of the Czarny Dunajec was drawn on the orthorectified aerial photo from 1964, and 31 river cross sections perpendicular to the centreline were delimited at 100-m intervals. The cross sections were subsequently replicated on orthorectified aerial photos from 1977, 1983, 1994 and the orthophoto from 2009, and measurements on the cross sections were done using the ‘Profile’ tool in Quantum GIS software. The lateral extent of the cross sections on all the river images was determined, accounting for lateral channel changes, and width of active channel, width of active river zone (including wooded islands within the river) and total width of low-flow channels in the cross sections were measured for all analysed dates. Low-flow channels represented parts of the active channel filled with water at base-flow conditions, during which the aerial photos were taken. On DEMs from particular dates, the following parameters were measured for each cross section: elevation of low-flow water surface, elevation of the lowest point of channel bed, mean elevation of active channel, mean elevation of active river zone and mean elevation of the belt of river migration in the years 1964–2009. With water surface in the mountain river frequently sloping obliquely to the cross sections, elevation of the water surface of low-flow channels was calculated as an average of the elevation of both low-flow channel margins, and in the case of a multi-thread channel, it was further averaged for all low-flow channels in a given cross section. Mean elevation of active channel, active river zone and the belt of river migration was calculated from measurements taken at 1-m intervals across the study cross sections.

While the procedure of extracting and editing of DEMs was sufficient to determine changes in the elevation of low-flow water surface and emerged parts of the study cross sections between the dates of given aerial photos, water depth in low-flow channels also had to be estimated in order to determine changes in mean elevation of channel bed. Water depth in the study cross sections was thus estimated on the archival aerial photos at 1-m intervals across each low-flow channel on the basis of photo-tone analysis and our knowledge of the river’s morphological conditions. Accuracy of these estimates was verified through repetition of this procedure for the orthophoto from 2016 and field measurements of water depth in the cross sections performed in 2017 at the same water stage at the Koniówka water-gauge station as when this photo was taken. These water depth measurements were also taken at 1-m intervals across low-flow channels in the cross sections.

Calculation of differences in vertical river position in the study reach between the dates of aerial photos required establishing limits of uncertainty related to elevation measurements in different parts of study cross sections. In the large sample of measurements taken in the cross sections (31 cross sections multiplied by a few tens to more than hundred measurements per cross section), random error component could be assumed to balance to a negligible cumulative value (cf. Schiefer and Gilbert 2007). Residual systematic errors observed for the network of check points after extraction and edition of a given DEM (Fig. 3b) resulted in an uncertainty component that was the same for all parameters characterizing vertical river position at a given date and ranged from ± 0.02 to ± 0.05 m for different dates of aerial photos. Further in the paper, we indicate procedure for correcting errors in the estimation of water depth in low-flow channels resulting from imperfect photo-tonal analysis of orthophotos. The analysed aerial photos/orthophotos were taken during periods of stable, clear weather, when river discharge was between 1.6 m3 s−1 and 3.3 m3 s−1. According to the rating curve for the Koniówka station, this range of discharges corresponds to 10 cm difference in water stage, and thus, the error in water depth estimation was assumed to be ± 0.05 m for all dates of aerial photos. This uncertainty limit thus referred to establishing the elevation of the lowest channel point at all dates, and with the proportion of low-flow channel width in the widths of active channel and river migration belt changing over the study period, it contributed between ± 0.02 and ± 0.025 m of uncertainty to the elevation of the former parameter and between ± 0.01 and ± 0.015 m to the elevation of the latter. Finally, the uncertainty related to imperfect editing of DEMs in forested areas was assumed to be ± 0.1 m. With an average proportion of forest in the study cross sections increasing over the study period from 11 to 45%, this factor added between ± 0.01 and ± 0.045 m to the overall uncertainty in establishing the elevation of river migration belt. As a result of the contribution of the three described factors, aggregated uncertainty differed between particular parameters characterizing vertical river position and dates of aerial photos, but for any date it did not exceed ± 0.1 m. The uncertainties in establishing changes in river elevation in a given period were the sum of uncertainties in establishing the river elevation at the beginning and the end of the period.

Data analysis

Some of analysed statistical samples have non-normal distribution, and thus, nonparametric statistical tests were used in the study. Statistical significance of the difference between estimated and measured water depth in low-flow channels was verified with a Wilcoxon signed rank test. Differences examined with this and other statistical tests were considered significant for p values < 0.05. A special approach had to be applied to the analysis of differences in river elevation in the study cross sections between the dates of aerial photos. First, as given cross sections were typified by consistent changes in vertical channel position (in terms of the direction of change) between successive dates of aerial photos, we limited the analysis to the years: 1964, 1983 and 2009. Choosing longer periods for the analysis improved the chance that changes in considered river characteristics between the beginning and the end of these periods will be meaningful and allowed for a better presentation of the results. Second, if the analysis was focused on the difference in the relative elevation of the study cross sections between 1964 and the later dates, all values in the initial sample would be zero and the lack of variance in this sample would preclude statistical analysis. Instead, a linear trend of river elevation in the study reach was estimated for 1964 on the basis of all cross sections, and river elevation in each cross section in relation to the elevation indicated by this trend was calculated for 1964, 1983 and 2009. Although longitudinal elevation profiles of rivers are typically concave upward on long distances (Wheeler 1979), in the 3-km-long reach of the Czarny Dunajec the linear model very well described a downstream decrease in river elevation (e.g. R2 = 0.996, p = 0.000001 for the linear trend of low-flow water surface). Linear trends estimated for 1964 were thus used as reference for calculation of the relative elevation of low-flow water surface, the lowest point of channel bed, active channel, active river zone and the belt of river migration at all analysed dates. Statistical significance of the differences in the so calculated relative river elevation between 1964, 1983 and 2009 was examined with a Friedman test, whereas the significance of differences between pairs of the years was verified with a Fischer’s least significance difference post hoc test. The same two tests were also applied to examine statistical significance of differences in river width between the dates of aerial photos.

Results

Estimated and measured water depth in low-flow channels

Estimation of water depth in low-flow channels shown on the orthophoto from 2016 and subsequent field measurements of water depth in the study cross sections indicated nearly the same cross-sectional variation of water depth (Wilcoxon test, p = 0.72), with the average coefficient of variation amounting to 0.549 for the estimated depth and to 0.541 for the measured depth (Fig. 5). However, mean values of the water depth determined on the basis of the orthophoto appeared to be underestimated by 0.12 m in comparison with the values derived from field measurements (Fig. 5), with the difference being highly significant statistically (p = 0.00001). These results suggested that the estimates of water depth derived from older aerial photos should be increased, while maintaining the same degree of water depth variation across the low-flow channels. To verify this, cross-sectional flow area in each study cross section was calculated as a product of low-flow channel width and mean estimated water depth and averages of the cross-sectional flow area in all cross sections were computed for each date of the archival aerial photos. These calculations indicated that the average cross-sectional flow area estimated from the 2009 aerial photo was considerably smaller (5.7 m2) than that estimated from older aerial photos (8.8 m2 in 1963, 9.3 m2 in 1983), and thus water depth estimates for 2009 were increased to yield a similar average value of cross-sectional flow area as for 1963 and 1983.

Fig. 5
figure 5

Mean and coefficient of cross-sectional variation of low-flow water depth in 31 cross sections of the Czarny Dunajec estimated from the orthophoto from 2016 and measured in the field in 2017. Statistical significance of the difference of the parameters between estimated and measured water depth, determined by a Wilcoxon test, is indicated

Changes in vertical river position in the study reach

Changes in vertical river position over the study period were determined for various parameters characterizing river position in the analysed cross sections (Fig. 6). Between 1964 and 2009, the elevation of low-flow water surface in the cross sections lowered by 1.57 ± 0.07 m on average and the change was statistically significant (Friedman test, p = 0.004) (Fig. 6a). In the years 1964–1983, the elevation of low-flow water surface significantly decreased by 0.95 ± 0.1 m (Fischer’s LSD test, p = 0.03), whereas its lowering by 0.62 ± 0.07 m between 1983 and 2009 was not statistically significant (Fig. 6a). The difference between the two periods is even more pronounced if annual rates of the lowering of low-flow water surface are compared—they equalled 5.0 cm/year in the first period and 2.4 cm/year in the second one. In 1964, deviation of the elevation of low-flow water surface in particular cross sections from the linear trend estimated for the study reach was relatively small, not exceeding ± 1.5 m (Fig. 6a). By 1983, it increased considerably as a result of channel incision and aggradation in the reach, with the former process prevailing over the latter. By 2009, the negative deviation from the linear trend increased to 3.6 m (Fig. 6a) indicating continued channel incision in the reach.

Fig. 6
figure 6

Changes of the relative elevation of low-flow water surface (a), the lowest channel point (b), active channel (c), and river migration belt (d) of the Czarny Dunajec River between 1964 and 2009. Boxplots show mean value (squares), median (the line inside boxes), quartiles (boxes) and extreme values (whiskers). Statistical significance of the difference of the parameters between the analysed dates, determined by a Friedman test, is indicated. Different letters denote statistically different samples indicated by a Fischer’s LSD post hoc test

Other parameters characterizing vertical river position in the study cross sections also changed significantly between 1964 and 2009. On average, the elevation of the lowest channel point in the cross sections lowered by 1.74 ± 0.17 m (Friedman test, p = 0.001) (Fig. 6b), that is more than that of low-flow water surface. Between 1964 and 1984, the elevation significantly decreased by 1.01 ± 0.2 m (Fischer’s LSD test, p = 0.01), but its further decrease by 0.73 ± 0.17 m in the years 1983–2009 was not statistically significant (Fig. 6b). In 1964, the deviation of the lowest channel point in the analysed cross sections from the linear trend estimated for the study reach did not exceed ± 1.2 m (Fig. 6b). By 2009, it expanded to ca. − 4.3 m as a result of channel incision in some parts of the reach and to 2.0 m as a consequence of bed aggradation in its other parts (Fig. 6b).

Average decrease in the elevation of the active river channel in the study cross sections between 1964 and 2009 equalled 1.54 ± 0.12 m and was statistically significant (Friedman test, p = 0.002) (Fig. 6c). This decrease reflected lowering of the active channel in the cross sections by 0.74 ± 0.14 m in the years 1964–1983 and by further 0.8 ± 0.12 m in the years 1983–2009 (Fig. 6c). In 1964, the deviation of the elevation of active river channel in particular cross sections from the linear trend determined on the basis of all cross sections ranged from − 1.50 to 0.76 m. By 1983, it increased to ca. ± 2.5 m, indicating a loss of vertical river stability in the study reach. By 2009, the negative deviation from the linear trend expanded to 3.8 m, whereas the positive deviation diminished to 1.65 m (Fig. 6c). Changes in the elevation of active river zone (including wooded islands) in the study reach were practically the same as those of active river channel.

Finally, the elevation of the belt of river migration in the years 1964–2009 lowered during this period by 1.03 ± 0.15 m on average, with the change being statistically significant (Friedman test, p = 0.016) (Fig. 6d). This total change reflected the lowering of the river migration belt by 0.42 ± 0.17 m in the years 1964–1983 and by 0.61 ± 0.17 m between 1983 and 2009 (Fig. 6d). In 1964, the elevation of the river migration belt in individual cross sections deviated by about ± 1 m from the values indicated by the linear trend determined for the study reach. By 2009, the positive deviation expanded to 1.83 m and the negative one to − 3.35 m (Fig. 6d).

In 1983 and 2009, median values of all the parameters charactering vertical river position in the study cross sections were lower than parameter means (Fig. 6). This indicates that channel incision occurred during the study period on a larger proportion of reach length than bed aggradation did. Indeed, such a spatial pattern of vertical channel changes is clearly visible in the longitudinal profile of the study reach, and Fig. 7 demonstrates such changes for the elevation of low-flow water surface. The figure shows that channel incision was recorded in the upper two-thirds of reach length, while bed degradation in the lower one-third. Moreover, differences in the temporal pattern of channel incision along the reach are apparent. Between 1.1 and 2.0 km of the reach length, the river incised mostly in the years 1964–1983, between 0.6 and 1.1 km similar amounts of incision were recorded in the periods 1964–1983 and 1983–2009, whereas in the uppermost 0.6 km of the reach most of incision occurred in the latter period.

Fig. 7
figure 7

a Extent of DEMs extracted from stereopairs of laterally overlapping aerial photos from 1964 and 1983 along the study reach of the Czarny Dunajec River. b Changes in the elevation of low-flow water surface along the study reach of the Czarny Dunajec between 1964 and 2009

Changes in river width

The measurements taken on the orthorectified aerial photos enabled us also to determine changes in river width between the consecutive dates of the photos. Even though this information was a by-product of the study focused on vertical channel changes, it appeared useful in explaining the latter changes, and thus, it is provided below. Width of the active channel in the study reach significantly decreased over the study period (Friedman test, p = 0.00001) (Fig. 8a). In 1964, it varied between 45 and 147 m, with the mean equal to 81.5 m. In 1983, it still highly varied between 25 and 127 m, but the mean significantly decreased to 51 m (Fischer’s LSD test, p = 0.00001). In 2009, the variation in active channel width decreased to the range 26–98 m, whereas the mean remained similar, amounting to 49 m (Fig. 8a). The width of active river zone in the study reach exceeded that of the active channel only slightly, by 1–8%, and the pattern of its temporal changes was similar.

Fig. 8
figure 8

Changes of active channel width (a) and low-flow channel width (b) of the Czarny Dunajec in the investigated river cross sections between 1964 and 2009. Boxplots show mean value (squares), median (the line inside boxes), quartiles (boxes) and extreme values (whiskers). Statistical significance of the difference of the parameters between the analysed dates, determined by a Friedman test, is indicated. Different letters denote statistically different samples indicated by a Fischer’s LSD post hoc test

Low-flow channel width in the study reach also decreased significantly between 1964 and 2009 (Friedman test, p = 0.04), but the timing of the change differed from that in active channel width (Fig. 8b). In 1964, the width varied between 15 and 61 m, with the mean amounting to 32.2 m. In 1983, neither the range of variation in low-flow channel width (14–58 m), nor its mean (29.5 m) differed significantly from the previous values. However, by 2009 mean low-flow channel width significantly decreased to 23 m (Fischer’s LSD test, p = 0.04) and the range of width variation diminished to 8.5–52 m (Fig. 8b).

Changes in sediment volume in the river migration belt

In the study reach, width of the belt of river migration in the years 1964–2009 varied between 50.9 and 157.5 m. Data about this width in successive river cross sections were combined with estimated change in mean elevation of the river migration belt in the cross sections to indicate changes in the volume of sediment removed (as a result of fluvial erosion or gravel mining) or deposited in the study reach during the years 1964–2009. This analysis was based on the assumption that changes determined for particular study cross sections are representative of 100-m-long river segments extending 50 m upstream and 50 m downstream of a given cross section.

Figure 9 presents results of the analysis for the whole period 1964–2009 and its earlier (1964–1983) and later (1983–2009) parts. During the years 1964–2009, an average of 11780 ± 1710 m3 of sediment was evacuated from a 100-m segment of the study reach, which corresponds to an annual loss of 256 ± 37 m3 of sediment per segment. Average values of the annual loss of sediment from a river segment during the periods 1964–1983 and 1983–2009 were similar and amounted to 267 ± 39 m3 and 248 ± 36 m3, respectively. However, while the upstream two-thirds of reach length experienced bed degradation and sediment loss, sediment deposition and bed aggradation occurred in the downstream one-third of the reach length (Fig. 9). In the upstream part, the average rate of annual loss of sediment equalled 453 ± 65.5 m3 per segment, whereas the average rate of annual sediment deposition in the downstream part amounted to 158 ± 23 m3 per segment. Both the disparity in these average rates and the considerably longer degraded part of the study reach explain the marked prevalence of sediment removal from the reach (Fig. 9).

Fig. 9
figure 9

Changes in sediment volume within the river migration belt in the study reach determined for successive 100-m-long river segments on the basis of DEMs extracted from archival aerial photos. Columns show a total change in sediment volume between 1964 and 2009, while solid and dashed lines indicate changes in the periods 1964–1983 and 1983–2009, respectively

The analysis also provided more detailed information about the temporal and spatial variation of changes in sediment volume within the study reach. In the downstream one-third of reach length, sediment deposition occurred in the earlier part of the study period, whereas in the later part no net change in sediment volume was recorded (Fig. 9). In the middle part of the reach, the main phase of sediment loss took place in the earlier part of the study period, whereas further upstream more sediment was removed in the later part (Fig. 9).

Discussion

The accuracy of photogrammetric measurements of the river channel depended on the quality of archival aerial photos (varying between the sets of photos from different years) and on field conditions. The photos were taken for different purposes: the flights from 1964 and 1983 for military applications, while those from 1977 and 1994 for civil ones. Storage conditions of the photos also varied, which was reflected in the occurrence of numerous scratches on the photos from 1964 and 1983.

The accuracy of aerotriangulation performed with the combined method (cf. Michalowska and Glowienka 2008) met the requirements for civil aerial photogrammetry established in Poland (GUGiK 1999) and was considered sufficient for the current study. For the aerotriangulation of the photos taken in 1964, 1977 and 1983, it was necessary to calculate parameters of camera distortion using ground control points, mainly located in test areas. For the photos taken in 1994, distortion parameters were known from a camera calibration certificate and no additional calculation was needed. Using a larger number of GCPs would allow for more accurate aerotriangulation, but as the RMS error of DEM extraction depends also on the flight altitude and pixel size, its value would be still at the range of 0.2–0.4 m for the analysed aerial photos (cf. Prokešová et al. 2010).

Automatic extraction and edition of DEMs for the valley of the Czarny Dunajec River required solving problems with uncertainty of DEM extraction in the areas covered by riparian forest and also in shadowed areas in case of insufficient radiometric resolution. Another cause of errors was the occurrence of specks or scratches on the images. All these factors caused that the extracted DEMs were the mosaic of properly and uncertainly measured areas. The latter usually coincided with riparian forest and had to be masked or edited. Similar problems caused by forest cover were met in the study of a landslide area in Slovakia (Prokešová et al. 2010). In such cases, the possibility for automatic measurement is limited and either edition of the DEM or visual measurement of the elevation is required (Dewitte et al. 2008).

Although the elevation of a given point on the DEMs may be indicated with some error, such errors are random rather than systematic (as inferred from their location and the values of residuals in bundle adjustment process), and thus, their influence on the accuracy of establishing vertical channel changes could be considerably reduced by the analysis of numerous river cross sections. Changes in the vertical position of low-flow water surface are considered to be the best indicator of vertical river changes modifying flow capacity of a channel that is channel incision or aggradation (Wyżga 1997; Wyżga et al. 2016c). The reach-average value of the lowering of low-flow water surface by 1.57 ± 0.07 m between 1964 and 2009 is considerably larger than the mean error of the DEM extraction for the study reach assessed by ground observations of the elevation of control points on the DEMs created for test areas (0.3–0.4 m), clearly indicating substantial incision of the Czarny Dunajec River during the analysed period.

The degree of uncertainty of the extracted DEMs is even less significant when compared with the amplitude of changes in the vertical position of low-flow water surface in the reach during the study period. The changes ranged from 3.61 m of channel incision in the middle and upper parts of the study reach to 1.34 m of channel aggradation in its lower part. In the middle part of the reach, channel incision exceeding 3 m since the mid-twentieth century was previously documented for selected river and floodplain cross sections, in which the elevation of paleochannels from different decades was compared with that of the contemporary river (Zawiejska and Wyżga 2010; Wyżga et al. 2016a, b, c). Here, the river not only dissected gravelly alluvium but also incised deeply into clayey bedrock and the strong lithological contrast between the gravels and the bedrock allows for easy comparison of the vertical position of paleochannels and the contemporary channel (Fig. 10). Channel aggradation in the lower part of the reach must have been facilitated by the delivery of bed material eroded in the incising upstream part of the reach, but its direct cause was the construction of a 1.8-m-high weir at the downstream end of the reach in the early 1960s (Fig. 11). Sediment deposition and channel aggradation in a short river stretch immediately upstream of the weir must have already occurred before the aerial photo was taken in 1964, and only minor changes in vertical channel position were recorded there during the study period. However, during that period considerable channel aggradation took place between 2.2 and 2.6 km of reach length (Fig. 7) as the river tended to adjust its longitudinal profile to the new erosion base created by the weir.

Fig. 10
figure 10

Deep channel incision in the middle part of the study reach caused dissection of non-resistant bedrock. The gravelly alluvium visible in the upper part of the undercut bank was deposited in the channel active in 1964

Fig. 11
figure 11

High weir constructed at the downstream end of the study reach in the early 1960s that has induced bed aggradation in the lower part of the reach

So far, only a few studies have reconstructed changes in vertical channel position on the basis of repeated surveys performed along a river course (e.g. Wyżga 1993; Piégay and Peiry 1997). The obvious reason for this is the scarcity of reliable longitudinal surveys of river channels from the past. The photogrammetric method offers a possibility of reconstructing vertical channel changes for different periods between the dates of archival aerial photos, which would provide information not only about the total amount of change, but also about its temporal and spatial dynamics. In the study reach of the Czarny Dunajec, the total amount of incision between 1964 and 2009 was similar along the upper two-thirds of the study reach, but the analysis of incision in the periods 1964–1983 and 1983–2009 clearly demonstrated the progressive upstream shift of incision. Backward erosion in the uppermost part of the study reach has long been retarded by the occurrence of resistant sandstone bedrock and the effective dissipation of flow energy in the multi-thread channel. However, artificial cutoff of one of the braids in 2005 concentrated river flow in the remaining single channel, and that induced its rapid incision (Wyżga et al. 2016a). This has accelerated backward erosion in the uppermost part of the reach, which now threatens the stability of a bridge located at the upper end of the reach (Fig. 12).

Fig. 12
figure 12

View of the uppermost part of the study reach where recent, deep channel incision has induced backward erosion threatening the bridge

Lowering of low-flow water surface reflects river incision increasing flow capacity of the channel, but channel deepening may also proceed without a lowering of low-flow water surface if it reflects the adjustment of channel geometry to a reduced sediment supply (Wyżga et al. 2016c). Channel incision and river metamorphosis may occur simultaneously—in that case the amount of channel deepening will be larger than that of the lowering of low-flow water surface. Such a situation was actually recorded in the study reach of the Czarny Dunajec, but the average lowering of the lowest point of channel bed exceeded that of low-flow water surface by only 10%. So small disparity between the two changes seems typical of rivers in the western part of the Polish Carpathians (cf. Wyżga et al. 2016c), which are mostly fed by non-cohesive, coarse-grained material (Wyżga et al. 2016b) and drain catchments with a limited degree of reforestation in the twentieth century. Much larger difference between the two changes was recorded in rivers draining the eastern part of the Polish Carpathians (Wyżga et al. 2016c), which flow through floodplains largely composed of cohesive sediments and experienced a considerable reduction in sediment supply during the twentieth century in response to reforestation of their catchments (Lach and Wyżga 2002).

Channel narrowing was only partly a spontaneous process of river adjustment to reduced water and sediment dynamics induced by in-channel gravel mining and land-use changes in the catchment (Zawiejska and Wyżga 2010; Wyżga et al. 2012). Channelization works were another important cause of the significant reduction in active channel width between 1964 and 1983. The works rapidly increased the proportion of reinforced channel banks in the study reach from 12% in 1964 to 34% in 1977, and since then, the proportion remained relatively stable (Hajdukiewicz and Wyżga 2019). The reduction in active channel width must have increased flow velocity and unit stream power at flood flows (cf. Wyżga 1993), hence stimulating further channel incision that has led to a reduction in low-flow channel width.

In highly energetic Polish Carpathian rivers, simultaneous occurrence or alternation of channel incision and lateral channel migration during the twentieth century has led to the formation of incised meander belts (Wyżga 2008; Wyżga et al. 2016c). Reconstruction of the lateral extent of channel migration on orthorectified archival aerial photos and of changes in the vertical position of channel migration belt with use of extracted DEMs is a valuable source of information for river management. First, it enables to determine a net change in sediment storage in a river reach, which is especially important in the case of reaches located upstream of dam reservoirs. In the study reach of the Czarny Dunajec, between 1964 and 2009 an annual rate of sediment loss from a 100-m segment equalled 256 ± 37 m3, and because large-scale gravel mining in the river was finished in the mid-1960s (Dudziak 1965), most of this volume must have been flushed out to downstream reaches. Second, it provides data for evaluation of the impact of the morphological change on hydraulic conditions in incised rivers. Rivers flowing in incised meander belts have very low potential for floodwater retention in floodplain areas because water is conveyed fast not only in the incised channel but also over low-lying floodplains (Wyżga 1999). Moreover, incised river reaches are typified by high shear forces exerted by flood flows on the channel boundary (Czech et al. 2016; Radecki-Pawlik et al. 2016), which is detrimental for the stability of the channels and in-channel hydraulic structures. So far, these hydraulic effects could be demonstrated through a comparison of river reaches with different morphological conditions, whereas the photogrammetric reconstruction of changes in three-dimensional river morphology offers a possibility of demonstrating temporal changes in hydraulic conditions caused by river incision.

Conclusions

In the last few decades many rivers worldwide underwent dramatic changes to their vertical position, with the direction and the rate of these changes depending on the type of impacts on the rivers and their catchments (Kondolf et al. 2002). Comprehensive knowledge about these changes along rivers is important for proper management decisions, but it can only occasionally be obtained through ground surveys (i.e. analysis of paleochannels) as rivers tend to rework their former channels. This is especially true for channelized rivers which can change their vertical position while remaining laterally stable.

Using the example of the Czarny Dunajec River in the Polish Carpathians, this study has demonstrated usefulness of DEMs extracted from archival aerial photos for reconstructing changes in vertical river position over the few past decades. Crucial for the final DEM accuracy was selecting numerous ground control points with stable position since the date of a given aerial photo, allowing to remove systematic errors of the DEMs. The analysis of sequential DEMs indicated considerable spatial variation in vertical channel change in the 3-km-long study reach—in its middle and upper parts, the river incised by up to 3.6 m, whereas up to 1.3 m of channel aggradation was recorded in its lower part. We were also able to detect temporal variation in channel incision along the study reach resulting in progressive upstream shift of bed degradation. Finally, combining the information about changes in vertical river position and width of the belt of river migration in the study period for successive river cross sections yielded data about the change in sediment volume along the river reach.

Twentieth-century incision of Polish Carpathian rivers recorded at water-gauge stations ranged from 0.5 to 3.8 m (Wyżga 2008; Wyżga et al. 2016c) and was larger than the accuracy of photogrammetric reconstruction of vertical channel changes. While vertical channel changes recorded at the stations may not be representative of longer river reaches, particularly in the upper course of the rivers, the direction, amount, spatial and temporal variation of recent vertical channel changes in river reaches outside water-gauge stations can be recognized through the analysis of DEMs extracted from archival aerial photos. As aerial photos from different dates since about the mid-twentieth century are available for rivers from many regions of the world, DEMs generated from these photos can be a useful tool in analysing vertical channel changes over that period.