Hydrological System of the Plitvice Lakes—Trends and Changes in Water Levels, Inflows, and Losses

Abstract

The Plitvice Lakes National Park is inscribed on UNESCO’s World Heritage List. The lake system is composed of 16 cascading lakes of different sizes separated by tufa barriers, which are the park’s key phenomenon. The lakes are characterized by highly diverse trends of the characteristic hydrological indicators—mean annual water levels, discharges, and tufa barrier growth. The analyses carried out in this paper identified that in the period before the early 1990s, Kozjak Lake had a trend of decreasing discharges, together with a trend of increasing water levels and growing tufa barriers. In contrast to this, in the period after 2001, a trend of increasing discharges was recorded, as well as a trend of decreasing water levels and decreasing tufa barriers. A potential cause of the barriers decreasing in size were the extremely high discharges during the last decade, which resulted in increased erosion of the tufa barriers. Losses of water due to the sinking from the lake system as well as the upper Korana course were confirmed, and it was identified that during the analyzed period the losses had not changed significantly. It was determined that the losses of water from Kozjak Lake occurred during low-water periods; however, they depended not only on the quantity of water flowing through the lakes but also on the hydrological conditions underground. The analyses carried out and the methodological procedures used in the analyzed area of the Plitvice Lakes are useful examples for the performance of analyses at similar lakes in karst formed by tufa deposition processes.

1. Introduction

Tufa precipitation and sedimentation processes in water systems cause the formation of tufa barriers with various morphological characteristics. The barriers formed in such environments are specific water systems, often recognized as very valuable phenomena and protected as unique natural systems. Many of these systems are inscribed on UNESCO’s World Heritage List [1], such as Pamukkala (Turkey), Jiuzhaigou (China), Yellowstone National Park (USA), Badad-e Surt (Iran), Band-e Amir (Afghanistan), Huamlong (China), the rivers of the Barkly karst (northern Australia), Hierve el Aqua (Mexico), Hanging Lake (USA), and Plitvice Lake (Croatia), which was the focus of the research presented in this paper. The formation of the barriers is largely connected with thermal springs and the chemical precipitation of tufa at locations with suitable physical conditions (temperature, water flow, etc.).

In some cases, the formation process of the tufa barriers is caused by very complex chemical and biological processes [2,3,4,5,6,7,8]. Many researchers have investigated various aspects of these aquatic systems, principally the conditions in which tufa precipitation occurs, as well as the speed of tufa precipitation, i.e., the formation dynamics of the tufa barriers. The focus has been, in general, on point observations of the growth at the chosen locations [9,10,11,12,13,14]. To the authors’ knowledge, estimation of the growth and erosion of the barriers based on the distributed hydrological estimations for the whole profile of the barriers has only been carried out at the Plitvice Lakes [15] and Krka River [16].

The formation and morphological changes of lakes in river systems are strongly influenced by the dynamics of tufa barriers’ growth. A relatively small area of the Dinaric karst in Croatia is rich in such systems. Well-known examples are the Zrmanja River, the Una River, the Mrežnica River, and, most significantly, the Krka River, which has large lake systems formed along it [17] and is recognized and protected as a national park in the same manner as the Plitvice Lakes.

The Plitvice Lakes represent an exceptionally valuable natural phenomenon lying in the central part of the Dinaric karst in Croatia. It is the oldest national park in Croatia and one of the oldest in Europe—it was established as a public institution in 1949. The process of tufa formation, which leads to the formation of tufa barriers and the creation of lakes, is the basic phenomenon of this national park and represents an outstanding universal value for which the Plitvice Lakes were inscribed on UNESCO’s World Heritage List in 1979 [1]. The national park covers an area of 296 km² [18], with the lakes accounting for about 1.9 km². As one of the most prominent researchers into the processes of the formation of the tufa barriers, who initiated the protection of the Plitvice Lakes in the form of a national park, Pevalek [19] stressed the specific quality of these processes: “Due to the action of tufa-forming plants, the Plitvice Lakes are an outstanding bio-dynamic system and it is only in a progressive state of that bio-dynamic system that the future of the Plitvice Lakes which we all wish for lies”. This paper deals with the hydrological aspects and the dynamics of the hydrological processes, with an emphasis on the unique characteristics of the Plitvice Lakes and the mutual conditioning of the hydrological relationships.

The lake system is composed of 16 cascading lakes of different sizes separated by the tufa barriers. The lakes lie roughly in the south–north direction over a distance of around 9 km, with an elevation difference of 162 m, between 637 m a.s.l. in the south and 475 m a.s.l. in the north [20]. At the end of the lake system, the karst river Korana starts its course, which is 134 km long [18]. It belongs to the Black Sea catchment area and drains into the Kupa River.

This lake system is characterized by dynamic processes and variable hydrological conditions, which have an impact on the change of the trends of some characteristic hydrological indicators as well as the status of the ecosystem associated with the water events. It’s dynamic nature is reflected primarily in the growth of the tufa barriers that raise and morphologically modify the lakes but also modify the hydrological characteristics of the water events in the impact area—i.e., water level fluctuations, inflows, and losses from the lake system [21]. For that reason, the immediate aim of the research carried out in this paper was to identify whether and how the hydrological interrelations related to the dynamics of the discharges, water levels, barrier growth, and water losses in the Plitvice Lakes study area changed over time, and then to provide an interpretation of the findings. The general objective was the selection of methodological procedures to analyze other, similar water systems in the karst regions where the processes of tufa barrier formation have a dominant role in the formation of lake-like systems.

2. Study Area

The Plitvice Lakes National Park lies in the mountainous karst region of inland Croatia, between Mala Kapela (1280 m a.s.l.), Lička Plješivica (1640 m a.s.l.), and Medveđak (884 m a.s.l.), at elevations ranging from 450 to 1280 m, latitudes of 44°44′34″ N to 44°57′48″ N, and longitudes of 15°27′32″ E to 15°42′23″ E. This is a carbonate bedrock region where the tufa extraction processes formed the lake system characterized by karst morphology and processes (Figure 1).

According to the data from the only climatological station in the Plitvice Lakes basin (44°53′ N, 15°37′ E) lying at an altitude of 579 m a.s.l., in the period from 1956 to 2019, the average total annual precipitation was 1615.5 mm, and the average annual air temperature was 8.8 °C. Very high variability is present, between both the precipitation levels for individual months and the levels of annual precipitation. The annual variations range between 1023.6 mm (2011) and 2213.9 mm (2014). The maximum and minimum monthly precipitation throughout the monitoring period were 447.2 mm (October 1974) and 3.0 mm (August 2000), respectively. The maximum and minimum mean monthly air temperatures were 21.7 °C (August 2003) and −9.3 °C (February 1956), respectively.

The lakes have a total volume of about 23 million m³ [23]. The lake system is divided into 12 Upper Lakes, which include the two largest ones, Prošće (area, 68 ha; volume, 7.67 million m³) and Kozjak (area, 82 ha; volume, 12.71 million m³), and 10 smaller lakes with a total area of 38.3 ha, and 4 Lower Lakes covering 6.4 ha in total [18,23]. The Upper Lakes mostly lie on the dolomite bedrock of poor permeability from the Triassic period, while the lakes downstream of Kozjak Lake and the Korana canyon lie on heavily permeable limestones from the Upper Cretaceous period [24,25,26] (Figure 1).

A characteristic of the Upper Lakes is that they are subject to the surface weathering and washout of the dolomites, which has an impact on the formation of wider valleys. On the other hand, the Lower Lakes cut like a canyon into the bedrock of permeable Upper Cretaceous limestone, the relative impermeability of which is provided by deep lacustrine deposits that fill up the space between the tufa barriers [21]. The Plitvice Lakes catchment belongs to the Black Sea catchment area, with its boundaries defined based on hydrogeological investigations and the control of the hydrological balance [20,22,27].

The largest quantities of water recharging the lakes come from the Crna Rijeka and Bijela Rijeka streams, which merge into the Matica River that flows into the most upstream lake—Prošće. In addition to the Matica River which, after flowing for around 600 m, discharges into Prošće Lake, a small tributary—the Sušanj stream—flows into the lake along its south-western shore. After flowing over the cascading tufa barriers of several smaller lakes (Ciginovac, Okrugljak, Batinovac, Veliko Jezero, Malo Jezero, Vir, Milino jezero, Galovac, and Gradinsko, Burgetići), the water from that lake system flows into Kozjak Lake, which has a significant permanent tributary—the Rječica stream. With the flow continuing, the water from Kozjak Lake flows through smaller lakes (Milanovac, Gavanovac, Kaluđerovac, and Novakovića Brod) and connects with the water from the left-shore tributary—the Plitvica stream—at the location of a semicircular depression called Sastavci (Figure 2a).

Sastavci is, at the same time, the place where the course of the Korana River formally starts. The early processes of tufa formation in that initial, canyon part of the course, where four so-called Korana waterfalls have also been formed, as well as in the more downstream section up to the so-called Korana bridge, are very intensive. In that part of the Korana course, karst forms—a network of sinkholes, caves, pits, and ponors—have developed on karstified limestone. All of this points to a highly developed secondary porosity and, consequently, intensive communication between the surface water and groundwater. It is precisely for this reason that, in spite of the marked processes of the formation of tufa barriers at certain waterfalls, new lakes have not yet formed, which would require making that part of the course impermeable through the processes of sedimentation [21].

When larger quantities of water flow in from the lake system, there is water flow present along the entire upper course of the Korana River. However, with decreasing inflows, the losses from the Korana channel increase, as a result of which certain sections of its upper course periodically run dry (Figure 2b). This regularly happens during prolonged dry summer periods. During such periods, there is no water in the river for a length of about 14.5 km, up until the permanent spring downstream of Vaganac village (Gavranića Vrelo), from where the Korana regains its permanent surface flow [20].

According to Hirc [28], who was the first to present the basic historical and geographical facts about the Plitvice Lakes, the lakes were first mentioned on a map of Croatia and Slavonia in the year 1673, with their first description provided by Demian in 1806. Franić [29] wrote the first detailed monograph about the Plitvice Lakes which focused mostly on the geographical aspects, mentioning that his research began in the year 1850, when a more detailed map was developed, the first lake depths were measured, and the first geological account was prepared. Koch [30,31,32] dealt with the hydrogeological and hydrographic conditions, and the first limnological investigations were completed by Gavazzi [33,34,35]. Pevalek [36,37,38,39] was the first to interpret the process of the formation of tufa and tufa barriers in the Plitvice Lakes.

In the 1950s, the first multidisciplinary scientific research project that included a detailed survey of the lakes was launched, continuous hydrological monitoring started, and the first hydrological analyses were carried out [40,41]. This period also saw the start or intensification of a number of other investigations documented in papers related to different aspects of the lake system, such as the geographic characteristics [42,43,44,45,46,47,48,49,50,51], the geological and hydrogeological characteristics [20,22,23,24,25,27,52,53,54,55,56,57], the properties of tufa and tufa-forming communities [5,58,59,60,61,62,63,64,65,66,67,68,69,70,71,72,73,74,75,76], the climate characteristics [77,78,79,80,81], the water and sediment quality, the dynamics of water exchange, and the biological and ecological issues of the lake system [82,83,84,85,86,87,88,89,90,91,92,93].

Since these were complex investigations, some of the papers mentioned above contained a hydrological component as well. Direct hydrological investigations, which are the very basis of this paper, were carried out by the authors of [15,21,94,95,96,97,98,99,100,101]. In the period 2016–2021, an interdisciplinary project under the title of “Hydrodynamic Modeling of Plitvice Lakes System” [102] was launched, involving the implementation of the recent hydrological investigations which are the subject of this paper and which were focused on studying the hydrological characteristics, trends and changes in water level fluctuations, inflows, and losses from the lake system during the period for which the measured hydrological data were available.

Information about the losses of water from the Plitvice Lakes has existed since the initial hydrological analyses and assessments. Three types of water losses have been identified: losses within the system of the barriers and lakes themselves, where some water is only diverted through a shorter underground flow through cavernous systems; real losses from the lakes towards remote points of discharge beyond the lakes themselves; and the most prominent type of water losses—from the course of the Plitvica tributary and the upper part of the Korana River course. For example, Franić [29] mentioned all the three types of losses, referring to the earlier papers by Gavazzi [34,35] and Hirc [28]. He also presented his own observations based on the accounts of the local population, which mentioned the sinking of water from Kozjak Lake. The losses from Kozjak Lake can be associated with the lake coming into contact on its right side with the permeable strata of carbonate structures, but they are still not fully direct due to the processes of deposition and precipitation of calcium carbonate in the lake region.

In order to detect the said losses, a series of discharge measurements were conducted periodically in the early 1980s at a number of profiles, mostly during dry periods when the losses were the most marked. This enabled the initial quantification of the hydrological correlations across the Plitvice Lakes, their tributaries, and the Korana River [103]. The problem of water losses in the analyzed basin area of the Plitvice Lakes and the upper Korana course was also dealt with by Beraković [104], Zwicker et al. [21], Biondić et al. [27], and Meaški [20], who linked them with the Kozjak fault.

The Plitvice Lakes system is characterized by spatial and temporal changes to the hydrological characteristics, both periodical, during intra-annual periods, and on a multiannual time scale. This was figuratively described in a paper by Petrik [41], who stated that a geodetic survey of the Plitvice Lakes in 1940 was carried out “due to the changes that take place in the lakes all the time without stopping”. The changes to the hydrological conditions affect the biogenic processes of the growth of the tufa barriers (Figure 3a,b) and their erosion (Figure 4a, which shows the concentrated overflow in part of the eroded tufa barrier), as well as the changes to a number of other ecological characteristics of that water system.

Almost all the natural processes that take place throughout the park are reflected precisely in the waterfalls and in the lakes. The fact that the lakes “grow” and that the tufa barriers are becoming taller and the waterfalls beneath them are becoming deeper has been observed already in a paper by Franić [29]. Petrik [41] analyzed the elevations of the lakes based on the geodetic surveys carried out by Omčikus in 1855 and Gavazi in 1912, and the results of the surveys from 1940 and 1951/1952 (the latter survey was related to a project that he had led in the 1950s). He concluded that even though the lake elevations were surveyed mostly in regard to their relative relations and different hydrological conditions, and the geodetic correlations were not too tight and were credible, “the barriers can grow significantly in multi-annual periods as we have illustrated and tried to prove on the example of Jezerce, but can also collapse due to breaking off, as seen beneath the entry barriers of Ciginovac and Okrugljak and in other locations”… It was also concluded that the differences in the results of elevation measurements (the results of geodetic surveys—author’s remark) shouldn’t therefore be considered a consequence of insufficient precision of surveying work, but rather a permanent characteristic of the lakes”. He stated that only in Gradinsko Lake were there three submerged barriers, two of which had been submerged long ago, while the third one had until a relatively recent date divided the lake into two parts—according to Franić [29], based on the results of Gavazi’s surveys, the elevation difference amounted to 8 cm. The results of the most recent survey were documented in a paper by Kapetanović et al. [44].

The uneven rate of growth is also reflected in the changing flows of water over the tufa barriers and the relative elevations of the lakes, an account of which is presented in the abovementioned paper by Petrik [41], as well as in this paper. The most striking example of the uneven growth is the 37 m high submerged barrier in Kozjak Lake (Figure 4b) that was submerged due to the faster growth of the downstream barrier, which is now the exit barrier from Kozjak Lake. According to Srdoč et al. [76], it was submerged some 430 years ago, before which time the flooded area of the present Kozjak Lake had been divided into two parts.

The rate of growth of the tufa periphyton on the Plitvice Lake barriers has for a long time been a subject of interest in a number of fields, but not in hydrology, at least not in such a way as for an attempt to be made to quantify the rate of growth of the tufa barriers based on the hydrological monitoring data. In other words, when analyzing the growth of the barriers or of the tufa itself in certain parts of the barriers, account has only been taken of the point data concerning the tufa growth at particular profiles or on artificial bedrock [105], and not of the growth of the barriers as a whole, which includes possible components of their erosion. As part of this paper, a hydrological assessment of the tufa growth rate was carried out using the example of an exit barrier from Kozjak Lake, based on an approach which was developed specifically in the Plitvice Lakes in a paper by Zwicker and Rubinić [15].

3. Materials and Methods

3.1. Available Materials

The systematic hydrological monitoring of the analyzed area started in the early 1950s, when the Prošće, Plitvički Ljeskovac, and Plitvice hydrological stations were established on Prošće Lake, the Matica River, and the Plitvica stream, respectively, in 1951 and the Kozjak Most hydrological station was established on Kozjak Lake in 1953. Later on, hydrological stations were gradually established on the other tributaries of the lakes. In the period 1991–1995, measurements at all the stations were interrupted due to the war. After 1995, the restoration of the hydrological stations started gradually, at first, mostly in terms of water level monitoring. As such, at the majority of the stations the discharges have only been monitored again since 2001.

Figure 5 gives a schematic representation of the position of the lakes, the hydrogeological cross-section of the Plitvice Lakes system, and the locations of the Plitvička Jezera climatological station as well as of several hydrological stations used in the paper.

Table 1. Basic information about the analyzed hydrological stations.

	Station	Stream/	Coordinates	Level “0”	Start of Water Gauge Operation/Restoration after Interruption	Analyzed Parameters	Analyzed Period 1
		Lake		(m a.s.l.)
1	Plitvički Ljeskovac	Matica	15°36′12″ N	636,412	August 1951/ November 1995	H, Q	1952–2020 1978–2020
			44°50′55″ E
2	Plitvički Ljeskovac	Sušanj	15°35′7″ N	639,091	September 1979/ January 2015	Q	1980–2020 2
			44°51′18″ E
3	Prošće	Prošće Lake	15°36′2″ N	636,487	August 1951/ September 1995	H	1952–2020
			44°51′54″ E
4	Plitvička Jezera	Rječica	15°36′52″ N	535,055	September 1979/ July 2002	Q	1980–2020
			44°52′22″ E
5	Kozjak Most	Kozjak Lake	15°36′29″ N	534,602 old 534,402 new	September 1953/ September 1995	H, Q	1953–2020 1953–2020
			44°53′35″ E
6	Plitvice	Plitvica	15°36′26″ N	556,132	July 1951/ January 1996	Q	1980–2020 3
			44°54′9″ E
7	Luketići	Korana	15°36′53″ N	403,967	September 1977/ November 1995	Q	1978–2020
			44°55′26″ E

¹ Starting from 1991, the monitoring was interrupted by the war, with the stations successively restored in the period from 1996 (Luketići–Korana) to 2002 (Plitvička Jezera–Rječica). ² Missing monthly data for the periods 2002–2014 and 2016–2020 were supplemented with the data from Kozjak Lake (correlation coef. 0.82). ³ Missing monthly data for 2014, 2016, 2017, and 2020 were supplemented with the data from the Korana–Luketići station (correlation coef. 0.85).

presents the basic information about the hydrological stations of the Croatian Meteorological and Hydrological Service (DHMZ) whose data was analyzed. It needs to be stressed that in addition to the interruption during the war, there is missing data in the data series for certain stations. In cases where it was possible to supplement this missing data with the data from the neighboring stations using correlation analyses, such an approach was used. Further on in the text, instead of the full names of the hydrological stations, the names of the streams/lakes on which the stations lie will be used. Likewise, due to the interruptions to the monitoring at all the analyzed hydrological stations which are documented in Table 1, the interruptions will not be expressed separately when discussing the individual analyses of the stations’ data which include the periods before and after the interruption of the monitoring.

3.2. Methodology

The time series of the hydrological data were analyzed by defining the basic characteristic statistical indicators (mean, standard deviation, minimum, and maximum) and the relation between the individual members of the same or compared series, starting from a homogeneity analysis, identification of the statistical significance of trends, and autocorrelation analyses.

The Wilcoxon test [106] was performed to assess the homogeneity of the input datasets and results. It is a nonparametric statistical test that compares two paired groups. The Wilcoxon test is performed by extracting an original series with n₁ members and a modified series with n₂ members from a data series of n = n₁ + n₂.

The measure to assess homogeneity was standard unit deviation U_0, which is defined using the equation:

$$U_0=\frac{S_0-E_{\left(S\right)}}{{\sigma }_S}$$

(1)

where:

S₀—modified series rank sum;

E_(s)—expected value of the modified series rank sum;

σ_s—standard deviation of the modified series rank sum.

The null hypothesis (that there was no significant change in the data series) is opposed by the alternative hypothesis (that there was a significant change, caused either naturally or in some other way). Taking into account a degree of confidence of α = 5% (which is common for hydrology) [106,107], the null hypothesis interval was adopted in accordance with the normal distribution:

$$-1.96\le U_0\le 1.96$$

(2)

The Mann–Kendall (MK) statistical test [108,109] has been widely used for detecting trends in hydrometeorological time series such as groundwater, water quality, streamflow, lake level, temperature, and precipitation [110]. The nonparametric Mann–Kendall estimate method was used to determine a monotonic increasing or decreasing trend for the variable of interest over time, with its statistical significance. The advantage of the nonparametric statistical tests over the parametric tests, such as the t-test, is that the nonparametric tests are more suitable for non-normally distributed, censored, and missing data, which are frequently encountered in hydrological time series [111,112].

The test is based on the correlation between the ranks and sequences of a time series. For a given time series {Xⁱ, i = 1, 2…, n}, the null hypothesis H₀ assumes that it is independently distributed, and the alternative hypothesis H₁ assumes that a monotonic trend exists in X [112]. The Mann–Kendall test statistic is expressed using the following equation:

$$S={\displaystyle \sum }_{k=1}^{n-1}{\displaystyle \sum }_{j=k+1}^{n}sgn(x_j-x_k)$$

(3)

where:

S—Mann–Kendall statistic, and

$$sng \left(x_i-x_k\right)=\{\begin{array}{c}1 if \left(x_i-x_k\right)>0\\ O if \left(x_i-x_k\right)=0\\ -1 if \left(x_i-x_k\right)<0\end{array}$$

(4)

The higher the absolute value of S, the stronger the proof that an increasing or decreasing trend exists. In case where there is a monotonic change in the time series, the expected value of S is 0, and the variance is calculated using the expression:

$$VAR \left(S\right)=\frac{1}{18}\left[n\left(n-1\right)\left(2n+5\right)-{\displaystyle \sum }_{p=1}^q{t}_p\left(t_p-1\right)\left(2t_p+5\right)\right]$$

(5)

where:

q—number of groups with the same values;

t_p—number of series members in group p with the same values.

In this paper, the Mann–Kendall tests are based on the calculation of Kendall’s tau as a measure of the association between two samples, which is itself based on the ranks within the samples. The computations assume that the observations are independent. Kendall’s tau values with a plus or a minus sign indicate a tendency of the analyzed parameter to decrease or increase, respectively. The null hypothesis was tested at a 95% confidence level for both water level and discharge data.

The codependence of the successive members of the same time series over a certain period was analyzed using an autocorrelation analysis. Based on the shape of the correlogram or the autocorrelation function (ACF), which is composed of a series of values of autocorrelation coefficients r_k, it was possible to identify periodicity in the occurrences of the analyzed hydrological series, as well as the system lag, i.e., its “memory”. According to Mangin [113], 0.2 was taken as a reference value, i.e., the lower threshold of the autocorrelation coefficient significance. Jukić [114] recommended that, since time series are non-stationary, time series as long as possible should be analyzed and that, for the computed ACF, it makes sense to analyze only the first, positive part of the function.

The losses of water from the lake system as well as from individual sections of the Korana River were analyzed using the water balance. In doing so, two types of discharge data series were used—those obtained from immediate discharge measurements in certain situations/on certain dates, and continuous data series concerning discharges at individual hydrological stations.

The dynamic of the tufa barrier growth/decline was analyzed according to the methodology described in the paper by Zwicker and Rubinić [15]. In doing so, use was made of information concerning the course of the mean annual water levels and discharges and their trends over a multiannual period, as well as the characteristics of the stage–discharge curve at the end of that period. The component of the average annual water level change was determined from the equation of the linear trend of the water level course, while the average annual change of the mean annual discharges was determined from the equation of the linear trend of the course of the mean annual discharges. Since the change in the discharge affects the water-level fluctuations, the overall change of the mean annual discharge value was determined from the analyzed trend of discharges. The stage–discharge curve from the reference hydrological station for the last year of the analyzed period was used to establish how much change the total change in the discharge brought to the water level in relation to the average water level; based on this, the component of the average annual water-level change caused by the discharge change was determined. The average barrier growth/erosion was determined by summing up the two components mentioned above, with a positive figure indicating barrier growth, and a negative figure indicating barrier erosion.

4. Results and Discussion

The lake system of the Plitvice Lakes is characterized by marked inflows in the section upstream of the highest-lying Prošće Lake, where the Matica River receives the water of the Bijela Rijeka and Crna Rijeka streams, which start their course from yield-rich springs. At the same time, the lake system also receives several concentrated inflows downstream of where the Korana starts its course (Figure 1). All the significant inflows were included in the hydrological monitoring that started back in the 1950s (Table 1). The run-off dynamics are also affected by the total volume of the lakes (more than 20 million m³), which makes the system almost lag behind the current precipitation and hydrological conditions, in relation to the intra-annual oscillations of the lake system and the possibility for a certain quantity of water to be retained in the active volume of particular lakes above the flow over the barriers.

In this context, the oscillations of water levels at the two biggest lakes (Prošće and Kozjak), as well as the discharges at the characteristic points of the lake system, were analyzed on both a multiannual and intra-annual time scale. Two components of the losses were analyzed—the losses from Kozjak Lake and from the Korana River. In this process, DHMZ data from the hydrological stations were used alongside data from several series of simultaneous discharge measurements which were carried out over two periods. In the period 1979–1987, 24 series of simultaneous measurements were carried out in different hydrological conditions in the Plitvice Lakes National Park. The measurement of losses stopped after 1987, only to be re-established again in 2015, with one simultaneous series per year. The losses of water from the Plitvica stream were not analyzed in this paper since detailed analyses of these have already been completed by Meaški et al. [56].

4.1. Precipitations, Air Temperatures, Water Levels and Discharges

The initial step in the analyses was to test the homogeneity of the characteristic climatic and hydrological indicators using the Wilcoxon non-parametric test [106]. The series of annual precipitation and the mean annual air temperatures (1956–2019) were analyzed, for which complete data series were available, partly supplemented with the data from the stations from the wider region. During the analyses, the periods up to and including the year 1987 and the period after that year were selected. It was identified that the analyzed precipitation series was homogeneous and that the air temperature series was non-homogenous. The reason behind the non-homogeneity was a highly marked trend of increasing air temperature during the analyzed period.

Analyses were also carried out of the homogeneity of the available data concerning the mean annual water levels at the stations in the two lakes, Kozjak (1953–1990 and 2001–2020) and Prošće (1952–1987 and 1997–2020), as well as the discharges for two data series—for the periods 1978–1990 and 2001–2020 at the Matica, Kozjak Lake, and Korana-Luketići stations. The obtained results showed that all the analyzed series of discharges were homogeneous and all those of water levels were non-homogeneous, the reason being the rising water levels in the lakes due to the growth of the tufa barriers. For this reason, an additional homogeneity test of the water levels within only the recent 20-year period was carried out, dividing the period into two identical parts (2001–2010 and 2011–2020). These results showed a period of homogeneous water levels.

Figure 6 presents the course of the mean, minimum, and maximum annual water levels at the Prošće (1952–2020) and Kozjak Lake stations (1953–2020) and the mean annual discharges for the Matica, Kozjak Lake, and Korana-Luketići stations (1978–2020), with marked trends for the mean annual values. The figure also presents the course trend of discharges in Kozjak Lake for an earlier period (1953–1990) and the total annual precipitation and the mean annual air temperatures at the Plitvička Jezera climatological station for the period 1956–2019. Even though it is methodologically questionable to define the trends for incomplete series of data about water levels and discharges, they are nevertheless presented here, since even to a small extent incomplete data series still show some general characteristics of the course of the analyzed data.

It is clear that in the period leading up to the rather long interruption of all monitoring activities during the 1990s because of the war, the mean water levels at Prošće Lake had a several times more marked increasing trend compared with Kozjak Lake, whereas after the re-establishment of monitoring, these trends were relatively close and slightly decreasing at both lakes. On the other hand, the mean annual discharges in the recent period, after the year 2000, showed a slightly increasing trend at all the stations. There were also clearly visible very high amplitudes of fluctuating values for the mean annual discharges due to extremely dry years such as 2011 or extremely wet years such as 2014. Furthermore, the last twenty years or so saw more frequent events with very high discharges, whose propagation through the lake system also affected the morphology of the tufa barriers in the form of increased erosion of the uppermost parts of the barriers.

Using the Mann–Kendall test [108,109], the trends of the mean annual water levels (divided into two periods: 1952–1987 and 1997–2020 in Prošće Lake; 1953–1990 and 2001–2020 in Kozjak Lake) and of the mean annual discharges for the 2001–2020 period were analyzed. The trends of the total annual precipitation and mean annual air temperatures (1956–2019) were also analyzed. Kendall’s tau is a non-parametric measure of correlation between the two ranked values. Kendall’s tau values with a plus or a minus sign indicate a tendency of the analyzed parameter to decrease or increase, respectively, and based on the resulting p-values, one can conclude whether or not a trend is statistically significant. The air temperature showed a statistically significant (p < 0.0001) increasing trend during the analyzed period (2.6 °C/100 years), whereas the precipitation showed no significant trend (p = 0.433) of change, with a very mild decrease of annual precipitation by app. 49 mm/100 years present.

The p values were lower than 5% only for the water levels, which meant that for the pre-1990 period, the trends for both stations were statistically significant, while for the later period, only the trend at the Prošće station was statistically significant. During the period ending in 1990, a marked linear trend (p < 0.0001) with an average annual water level increase of 1.12 cm was identified at the Prošće station. At the Kozjak Lake station, the trend was slightly milder but still significant (p = 0.003), with an average annual water level increase of 0.37 cm. However, after the year 1997, the situation was reversed, with a visible negative linear trend with an average annual water level decrease of 0.29 cm at the Prošće station (p = 0.049) and 0.19 cm at the Kozjak Lake station (p = 0.398).

From the analysis of the discharges, the p-values exceeded 5%, which meant that even though trends existed, they were not statistically significant. It was noticeable that the trends were on the rise at all the three stations, with an average annual discharge increase of 0.019 m³/s at the Matica station (p = 0.770), 0.036 m³/s at the Kozjak Lake station (p = 0.581), and 0.011 m³/s at the Luketići-Korana station (p = 1.0).

The mean monthly and annual water levels in both lakes were also analyzed (

Table 2. Mean monthly and annual water levels (cm) at the Prošće and Kozjak Lake stations (2001–2020) (Note: Avg—average; SD—standard deviation; Cv—coefficient of variation; Max—maximum; Min—minimum; MAX—absolute maximum; MIN—absolute minimum).

	January	February	March	April	May	June	July	August	Septemebr	October	November	December	Year
Prošće Lake
Avg	68	69	73	74	70	63	57	53	55	61	67	70	65
SD	7	7	7	9	8	6	5	7	9	10	9	9	5
Cv	0.10	0.10	0.09	0.12	0.12	0.10	0.10	0.13	0.17	0.16	0.13	0.13	0.07
Max	79	80	82	87	83	75	64	71	87	81	80	81	75
Min	57	54	60	59	54	51	45	44	46	47	48	52	55
MAX	95	95	100	101	111	91	89	87	103	110	98	97	111
MIN	49	53	55	55	52	47	42	42	42	44	46	47	42
Kozjak Lake
Avg	101	104	110	111	105	96	88	84	85	89	97	103	98
Stdev	11	12	12	15	12	10	8	9	13	14	12	14	8
Cv	0.11	0.11	0.11	0.14	0.11	0.10	0.09	0.11	0.15	0.15	0.13	0.14	0.08
Max	122	126	127	137	128	118	102	108	131	119	115	123	112
Min	75	72	83	85	84	80	76	70	69	68	66	72	81
MAX	156	150	162	161	187	138	129	129	162	167	145	154	187
MIN	71	72	75	77	81	76	73	68	67	64	65	65	64

). Due to the non-homogeneity of the multiannual series of the water level data, Table 2 presents the results of that analysis for only the last 20-year period of their joint operation (2001–2020), which were demonstrated as homogeneous by the test results.

In Prošće Lake, in the analyzed recent 20-year period of observation (2001–2020), the mean annual water level was 65 cm, ranging between 55 cm and 75 cm. During that period, the absolute maximum and minimum levels were 111 cm and 42 cm, respectively. The amplitudes of the annual water level oscillations ranged between 31 and 61 cm, with an average of 46 cm. In that same analyzed period, the mean annual water level in Kozjak Lake was 98 cm, ranging between 81 and 112 cm. The absolute maximum and minimum levels were 187 cm and 64 cm, respectively. The amplitudes of the annual water level oscillations ranged between 43 and 109 cm, with an average of 69 cm. As for the occurrences of extreme values during the year, it was identified that they mostly occurred during the same months in both lakes.

Table 3. Mean monthly and annual discharges (m³/s) at the Matica, Kozjak Lake, and Korana-Luketići stations (1978–2020) (Note: Avg—average; SD—standard deviation; Cv—coefficient of variation; Max—maximum; Min—minimum).

	January	February	March	April	May	June	July	August	September	October	November	December	Year
Matica
Avg	2.26	2.47	3.34	3.78	2.83	1.82	1.29	1.05	1.26	1.69	2.20	2.64	2.22
SD	0.946	1.11	1.61	1.78	1.36	0.678	0.425	0.467	1.21	1.16	1.07	1.30	0.61
Cv	0.42	0.45	0.48	0.47	0.48	0.37	0.33	0.44	0.97	0.69	0.48	0.49	0.28
Max	4.68	5.01	7.34	8.19	5.89	3.45	2.57	2.85	7.66	5.44	4.33	6.36	4.22
Min	1.01	0.796	0.859	1.25	0.967	0.824	0.601	0.577	0.607	0.648	0.689	0.731	1.25
Kozjak Lake
Avg	2.85	3.21	4.14	4.89	3.96	2.63	1.78	1.43	1.56	1.89	2.50	3.19	2.84
SD	1.19	1.43	2.00	2.11	1.78	0.939	0.588	0.691	1.40	1.42	1.30	1.48	0.813
Cv	0.42	0.45	0.48	0.43	0.45	0.36	0.33	0.48	0.89	0.75	0.52	0.46	0.29
Max	5.35	6.61	8.19	9.06	8.29	4.74	3.30	4.04	8.79	6.61	5.81	6.95	5.58
Min	1.00	0.960	1.35	1.91	1.27	0.938	0.775	0.755	0.678	0.674	0.647	0.905	1.42
Korana-Luketići
Avg	2.61	2.99	4.34	5.16	3.84	2.14	1.07	0.729	1.15	1.54	2.19	3.09	2.57
SD	1.38	1.71	2.25	2.45	2.29	1.17	0.624	0.726	1.97	1.79	1.44	1.87	0.89
Cv	0.53	0.57	0.52	0.48	0.60	0.55	0.58	1.00	1.71	1.16	0.66	0.61	0.34
Max	5.93	7.55	9.01	10	9.86	4.65	2.33	3.31	11.3	7.25	5.10	8.41	5.27
Min	0.703	0.479	0.641	1.03	0.590	0.580	0.146	0.032	0.011	0.114	0.020	0.230	1.050

presents the mean monthly and annual discharges at the Matica, Kozjak Lake, and Korana-Luketići stations for the joint period of observation 1978–2020. The mean annual discharge at the Matica station at the start of the inflow into the lake system was 2.22 m³/s; at the central station of the lake system, Kozjak Lake, it was 2.84 m³/s; and at the most downstream station of the lake system, Korana–Luketići, it was 2.57 m³/s. It was clear that the distribution of discharges at the section between the Matica and Kozjak Lake stations was standard—for every month, the more downstream station had higher average discharges due to the contribution of the tributaries from the intermediate basin, with the intermediate basin contributing to the balance with 0.62 m³/s on the annual level. As opposed to this, the value of the mean annual discharge at the Korana-Luketići station was 0.27 m³/s lower than that at the Kozjak Lake station lying more upstream, despite the fact that, in addition to the immediate intermediate basin, water from the significant tributary, the Plitvica stream, also flows into that initial section of the Korana course. The mean annual discharge of 0.667 m³/s was measured at the Plitvica Lake over a slightly shorter available period of discharge data (1980–2020). Hence, the difference in discharge was even more marked, ranging between 0.9 and 1 m³/s. This suggested that there were considerable losses in the upper Korana course, caused by the geological composition of the bedrock. The difference in the average monthly discharges between Korana-Luketići and Kozjak Lake was positive (higher inflows at the more downstream station) only during March and April, when the mean monthly inflows were, on average, the highest. All the other calendar months saw lower discharges in the more downstream profile due to the abovementioned losses.

In Kozjak Lake, 78% of the inflow was accounted for by the water flowing into the Plitvice Lakes system through the Matica watercourse and 16% through the Rječica watercourse, whereas the total contribution of the Sušanj watercourse and the intermediate basin accounted for app. 6%. It is possible that this inflow was slightly higher, since in Kozjak Lake there were losses in certain parts of the year, as well as a water supply intake structure (55 L/s during the summer months and 30 L/s during the other months), and their balance share was not contained in the registered discharges (overflow water) in Kozjak Lake, since there was no appropriate monitoring of these losses. Only around 73% of the total average discharges from Kozjak Lake and Plitvica Lake flowed into the Korana-Luketići profile, and 27% of the water was lost through infiltration into the channel.

The autocorrelation characteristics of the course of daily discharges at the studied stations were analyzed as well (Figure 7). The diagram displays marked seasonality in the course of daily discharges at all the three stations, as well as their mutually well-balanced course. If, as according to Mangin [113], 0.2 was taken as the reference value, i.e., the lower threshold of the autocorrelation coefficient significance, the Kozjak Lake station lagged the most. At that station, the ACF had a value higher than 0.2 in its initial falling limb over 71 days, whereas at the Korana-Luketići and Matica stations the duration of codependence was slightly shorter, amounting to 55 and 56 days, respectively. This indicated the impact of the mass of water retained in Kozjak Lake on the regulation of discharge from that lake.

4.2. Losses from Kozjak Lake

The losses from Kozjak Lake were not as easily discernible as the losses at the most upstream section of the Korana course that make the Korana run dry for several dozen days a year. These losses are mutually dependent, since on the one hand, the water losses are affected by the hydrological conditions of the surface water and the conditions in the underground, and on the other hand, the losses of water into the underground affect the hydrological conditions of the surface water.

Still, the losses from Kozjak Lake did exist, and were even mentioned by Franić [29], who referred to the measurements carried out by the engineer Setinski in 1909. During the construction of a water pipeline from the Rječica spring, Setinski also measured the discharges at a number of profiles along the lake system of the Plitvice Lakes during the dry period. By performing these measurements, he, among other things, identified a discharge of 1162 L/s from Prošće Lake and 1043 L/s from Kozjak Lake. Although there were certainly tributaries that were unmeasured on this occasion, there was a difference of around 119 L/s between the lakes, which was a definite confirmation of the losses of water from the lake itself. Similar results were obtained by Setinski’s measurements at the more downstream section of the course through the Lower Lakes and in the most upstream part of the Korana course. Petrik [41] analyzed these results, questioning the possibility of the water sinking within the lake system of the Plitvice Lakes. He concluded that the only signs of the underground flow of water in the Plitvice Lakes were found in the tufa, and that this was water circulating within the tufa barriers. Meaški [20] linked the possible losses of water from the lake system of the Plitvice Lakes with the run-off of water shifting from the area of poorly permeable dolomites to the well-permeable area on the north-eastern side of the Kozjak fault (Figure 1).

In this paper, the balance of the losses of water from Kozjak Lake was determined using two types of data—analysis of the inflows–losses correlation based on the processed data concerning the mean daily and mean monthly discharges from DHMZ’s Kozjak Lake, Matica, Rječica, and Sušanj hydrological stations, and based on the simultaneous results of DHMZ’s discharge measurements at the Kozjak Lake, Prošće, and Rječica profiles.

If the losses of water from Kozjak Lake were analyzed based on a simplified approach—the results of simultaneous discharge measurements at the profiles (Figure 8)—it would be clear that the losses from the lake on average occurred when the discharge in Kozjak Lake was around 1.6 m³/s or lower.

The assessment of losses was also made by analyzing the inflows–losses correlations based on the data from DHMZ’s Kozjak Lake, Matica, Rječica, and Sušanj hydrological stations for the joint period 1980–2020. It needs to be stressed that due to the missing data, these analyses do not include all the water-balance parameters, i.e., precipitation, evaporation, inflows from the intermediary basin, and abstraction. These are the balance elements that have plus/minus signs in opposition to each other (precipitation–evaporation, inflows from the intermediary basin–abstractions), for which the preliminary assessments of their balance shares identified that they do not significantly affect the overall results of the assessment of losses. For this reason, the losses were analyzed only on the basis of the conducted measurements of inflows into Kozjak Lake from the Matica, Sušanj, and Rječica hydrological stations, and of the outflow from the lake at the Kozjak Most hydrological station.

Figure 9 presents the mean monthly discharges from the Kozjak Lake hydrological station and the discharge difference (= Kozjak Lake − (Matica + Sušanj + Rječica)) for two periods, 1980–1990 and 2002–2020.

It was clear that during the first period, losses occurred to a greater extent, with a negative discharge difference for as long as 6 months. In the period 2002–2020, the discharge difference was negative during only three months—September, October, and November. In both periods, the heaviest losses were recorded after the summer, i.e., during October, with an average of 0.275 m³/s for the 1980–1990 period (13% in relation to the inflows) and an average of 0.146 m³/s for the 2002–2020 period (7% in relation to the inflows), or 0.193 m³/s overall (9% in relation to the inflows).

If the monthly losses during all the years in the analyzed period were studied, it was only in 2005 (which was particularly wet during the dry part of the year) that losses were not recorded throughout the year. Such distribution of the water discharge difference, with the highest negative differences during October and November, indicated that the losses were affected not only by the hydrological conditions of the surface water (water level and inflows into Kozjak Lake), but also by the conditions underground. That is, the months with the highest average losses did not have the lowest discharges. One can therefore assume that the higher losses in early autumn were the consequence of the slower reaction of the rising groundwater levels surrounding Kozjak Lake in relation to the inflows into the lake system, and consequently of the higher hydraulic gradients of the underground run-off from the lake. Regretfully, for the time being, there are no results of groundwater observations which might be used for a more detailed analysis of the problem of water losses from Kozjak Lake, but it can be expected that this will be one of the most important reasons for the establishment of such observations.

4.3. Losses in the Upper Korana Course

The losses of water from the Korana channel in its upper course were clear. The losses were of such a nature that during low-water periods, the upper part of the Korana course frequently, even completely, ran dry downstream of Sastavci, where constant flow was ensured by the water flowing over from the Plitvice Lakes and the inflow from the Plitvica stream. The occurrence of losses was also mentioned by Hirtz [28], Gavazzi [34,35], Franić [29]; Petrik [41], DHMZ [103], Beraković [104], and Zwicker et al. [21] offered their initial quantifications for the losses. In the paper by Biondić et al. [53], the focus was on identifying an underground link between the water in the Korana catchment and the Klokot spring in neighboring Bosnia and Herzegovina. In order to accomplish this, the result of the Rastovača pit dye tracing conducted in 2005 was analyzed, confirming the existence of such a link. This also confirmed that the waters that sink from the upper Korana course run off towards the neighboring Una watercourse, where they appear in the form of its left-bank springs in the wider area of the Klokot spring zone near Bihać.

Just as when analyzing the losses from Kozjak lake (Section 4.2), the losses of water from the upper part of the Korana course were also analyzed using two types of data—analysis based on a series of simultaneous discharge measurements by DHMZ during low water periods, and analysis of the inflows–losses correlations based on the data from DHMZ’s Kozjak Lake, Plitvica, and Korana-Luketići hydrological stations. In the period from 1979 to 1987, simultaneous measurements were carried out in different hydrological conditions at six profiles along the Korana River—at the Kozjak Lake and Plitvica stations, which represent the water balance at the start of the Korana course, Lajtmanov Slap, Luketići, Korana Most, Magdića Mlin, and Čatrnja. Simultaneous discharge measurements started to be taken again in 2015 on a yearly basis, but from a slightly increased number of profiles. However, for comparison purposes, when processing data, this paper analyzed only the results of the discharge measurements at the profiles which were included in periods: the start of the Korana course (Kozjak Lake and Plitvica), Lajtmanov Slap, Luketići, and Korana Most.

Figure 10 presents the position of the selected profiles and the comparison of the mean discharges from all the discharge measurements for the period 1978–1987, with the losses between each profile for the two abovementioned periods, 1978–1987 and 2015–2020. It is clear that the losses at the Korana River were very marked. They occurred between every discharge measurement profile in both periods. The average discharge decrease, between the Kozjak Lake and Plitvica profile and the Korana Most profile, was 1.12 m³/s (63% of the losses, in relation to the discharges measured at the initial Kozjak Lake and Plitvica profile) in the period 1979–1987, and 1.32 m³/s (as much as 89% of the losses) in the period 2015–2020. Up to the Luketići profile, the average discharge decrease was 0.93 m³/s (52% of the losses) in the period 1979–1987, and 0.94 m³/s (63% of the losses) in the period 2015–2020, with a total average of 0.93 m³/s (55 %).

It is obvious that during the two relevant periods the losses had not changed, i.e., that the characteristics of the sinking zones in the analyzed section of impact of the Korana course had not changed.

The relation between the losses and discharges measured at the Luketići profile for the 1978–1987 and 2015–2020 periods during the attainment of a series of discharge measurements is presented in Figure 11. The figure makes it clear that the relation between the losses and discharges was very strong (coefficient of determination R² = 0.96), confirming that the results from the recent period matched the results from the old period very well, which meant that there were no changes in the characteristics of sinking in the analyzed section of the course. It is also clear that the lower the discharge, the higher the rate of losses. The maximum loss of 100% was measured in October 1983, which meant that at the Luketići profile, there was no water at all. The minimum loss of 20% was measured in April 1982, which was at the same time as the measurement for the highest discharge (6.94 m³/s at Kozjak Lake and Plitvica).

Figure 12 presents the total number of days and the number of consecutive days for each analyzed year from 1978 to 2020 when the Korana River ran dry at the Luketići hydrological station. The year 2003 had the highest total number of days with the Korana River running dry, 57 days, and the highest number of consecutive days, 36. For the same period and the same station, Figure 13 presents the total number of days per month when the Korana River ran dry. One can see that September, with a total of 128 days or an average of 11%, was the month with the highest number of running-dry events during the analyzed period.

4.4. Tufa Barrier Decrease Estimation

The formation of the Plitvice Lakes was affected by the formation of tufa barriers which dammed the original river valley, leading to the formation of the lakes. The growth of the tufa barriers is a biodynamic process that has taken place for thousands of years in conditions of undisturbed natural relations. Each barrier has its own dynamics—some grow at a faster rate, others at a slower rate. In this process, the downstream barriers sometimes rise above the upstream ones which then become submerged, with two lakes merging into one (Figure 4b). There are also instances of the barriers decreasing in height because of their erosion and breaking off due to falling trees or ice melt after the winter period, which has an impact on the change in the flow and dynamics of the water. The growth of the tufa barriers raises the water level and changes the volume of water in the lake. In addition to the barriers, calcium carbonate extracts in the lakes themselves form deposits on the bottom, creating fine lacustrine mud. In this way, the lake is slowly filled up and the morphology of the bottom is modified [100].

The trends of individual hydrological values are an important indicator of the hydrological and morphological processes in water systems. Because of their connection, changes in the trends of hydrological values also indicate changes in the morphological processes. As was mentioned in the introduction, the tufa barriers of the Plitvice Lakes have a general growth tendency, both during the historical period and in recent times. In a paper by Srdoč et al. [76], it was determined that the average rate of increase of the lake level was as high as 1.35 cm per year, or around 17 times faster than the rate at which the lake bottom became covered with mud. Similar estimates about the barrier growth rate (on average, 13 mm per year) were also presented in a paper by Božićević et al. [43]. A much wider range for the possible barrier growth, on average 1–3 cm per year, was given by Petrik [41], based on the estimated differences in the geodetic surveys of the lakes’ positions, which had been carried out a century earlier by Omčikus (in 1855), who noted that the elevation position of Kozjak Lake changed the least during that time.

The first quantification of the growth of tufa barriers in Kozjak Lake using a hydrological approach was presented in a paper by Zwicker and Rubinić [15] for the period 1952–1990, i.e., for the period up until the temporary interruption of monitoring. For the analyzed 39-year period, it was determined that the average annual growth of the Kozjak Lake barrier was 0.56 cm/year [15]. This indicator of the average growth included the components of both barrier growth and of incidental barrier erosion, and did not represent point data but referred to the entire overflow profile across the tufa barrier.

For the period after the re-establishment of monitoring in 2001, significant changes were observed in the trends of water levels and discharges (Figure 6) [95,96,115]. The differences in the abovementioned papers were also more marked partly due to an error at level “0” of the restored Kozjak Lake station. Even after the correction of the elevation changes based on the more recent geodetic survey, the differences in the courses of the mean annual water levels of Prošće and Kozjak Lake and the discharges in Kozjak Lake were clear (Figure 6)—the mean annual water levels in both lakes in the recent period, after the monitoring was re-established, showed a change to a decreasing trend. The mean annual discharges also showed a change of trend: after having a decreasing trend before the monitoring was interrupted, they now have an increasing trend.

In the same way as in the paper by Zwicker and Rubinić [15], the hydrological conditions and trends of the mean annual water levels and discharges of Kozjak Lake were analyzed in the period after 2001. It was identified that the mean annual water levels during the analyzed 20-year period (2001–2020) had a summary decrease of 4 cm, and the mean annual discharges had an increase of 0.722 m³/s. According to the stage–discharge curve defined by DHMZ for the year 2020, this led to a summary change of 7 cm. The total change was 11 cm, which gives an average annual decrease in the tufa barriers of 0.55 cm/year, which is practically identical (aside from having the opposite sign) to the figure calculated for the average growth of the Kozjak Lake tufa barrier for the period up until 1990.

The alternating periods of the growth and erosion of the tufa barriers are natural processes, and the process of barrier decline rather than growth has dominated in recent times. The likely reason for this is the erosion of barriers due to the observed trend of increasing discharges. In other words, high-water levels are one of the dominant triggers for barrier erosion and destruction, as was the case in April 2018. On 17 March 2018, there was a high water wave with a peak discharge of around 19.3 m³/s at the Kozjak Lake station, with such high water levels last recorded as far back as the year 1976, when a discharge of 25.1 m³/s was recorded. After only 14 days (1 April 2018), there was another high water wave with a peak discharge of 18.9 m³/s at the Kozjak Lake station. As the result of two such extremely high water waves over a short period, a barrier at Milino Lake broke off, followed by strong erosion (Figure 14), after which the conditions of run-off and discharge in the lake changed, with the lake itself drastically reducing its normal water surface.

Figure 15 presents in detail the course of the mean daily discharges at the Matica and Kozjak Lake stations for the period from 2018 to 2019. The chart marks the barrier breaking off on 9 April 2018. It also shows that only a year later, a much higher water wave with a peak discharge of 29.8 m³/s at the Kozjak Lake station was recorded, which was the highest water wave recorded so far during the monitoring period. Frequent high-water events have caused the abovementioned changes in the trends of the annual values of water levels and discharges in the other parts of the lake system as well. The growth and erosion of the tufa barriers and the related changes in the dynamics of the fluctuation of water levels in the Plitvice Lakes are a natural process, but sudden changes in the lakes’ morphology are nevertheless undesirable in an area of the national park that has a large number of visitors.

5. Conclusions

The analyses made in this paper showed the marked complexity of the functioning of the hydrological system of the Plitvice Lakes, a valuable site on UNESCO’s World Heritage List. The study identified that the regime of fluctuating water levels in the lakes and the flow of water through the lakes was also affected by the morphological changes in the geometry of the lakes themselves due to the variable growth of the tufa barriers that form the lakes. Since the monitoring of hydrological conditions in the analyzed area was fully interrupted in the 1990s due to the war, this paper compared the annual trends of data on the fluctuations of the characteristic annual water levels and discharges before and after said interruption. In addition, the courses of the annual precipitation and the mean annual air temperatures were also analyzed.

This water system is characterized by a mean annual precipitation of 1615.5 mm, a mean annual air temperature of 8.8 °C, and a mean annual discharge of 3.24 m³/s. During the analyzed period of 1956–2019, the annual precipitation showed no significant trend of change (a very mild decrease in the annual precipitation by app. 49 mm/100 years was recorded), while the course of the mean annual air temperatures showed a very statistically significant upward trend of 2.6 °C/100 years, which can be attributed to the manifestations of climate-change impacts.

It was determined that the data concerning the annual trends of the water levels in the two biggest lakes (Prošće and Kozjak) from the periods before and after the monitoring was interrupted were not homogeneous, whereas the data concerning the trends of the analyzed discharges were homogenous. The reason for this was the statistically significant trends of water levels in the lakes. Namely, in the period up until the early 1990s, Kozjak Lake had a trend of decreasing discharges and increasing water levels (0.369 cm/year). In contrast, the period after 2001 recorded a trend of increasing discharges and decreasing mean annual water levels (0.189 cm/year). It is clear that the trends of water levels and discharges changed their course after the monitoring was interrupted. All of this affected the growth of the tufa barriers, whose rate of growth recently became negative, after an earlier average annual increase of 0.56 cm/year. There was also an average annual decrease in the water level in Kozjak Lake by 0.55 cm/year. The reason behind this could be the increasing frequency of extremely high-water events, such as the one which broke off part of the barrier in Milino Lake in 2018. In addition to having an impact on such extreme changes in the lakes’ morphology, the high-water events also affect the accelerated erosion of the upper parts of the barriers over which the water flows.

The paper also analyzed the losses of water both from Kozjak Lake and from the upper part of the Korana course. Even though the losses of water from Kozjak Lake had been identified more than a century ago [29], doubts about their occurrence have been present for a long time. According to Petrik [41], these were not real losses but rather the circulation of water within the tufa barriers. Meaški [20] linked the losses of water from the lake system of the Plitvice Lakes with the water run-off shifting from the area of poorly permeable dolomites to the well-permeable area on the north-eastern side of the Kozjak fault. This paper confirmed the existence of the losses and estimated that, on average, they occurred when the discharge in Kozjak Lake was around 1.6 m³/s and lower. Their occurrence was most often associated with the period at the very end of the dry season (September–November), when the autumn surface inflows already occur as the consequence of autumn rains. The average amount of losses from Kozjak Lake of 0.193 m³/s constituted 9% of the total inflow. Losses were recorded in both analyzed periods: for the otherwise drier 1980–1990 period, an average of 0.275 m³/s (13% in relation to the inflows), and an average of 0.146 m³/s for the 2002–2020 period (7% in relation to the inflows).

The identified intra-annual distribution of losses indicated that the losses were affected not only by the hydrological conditions of the surface water (water levels and inflows into Kozjak Lake), but also by the conditions underground. Namely, the months with the highest average losses did not have the lowest discharges. One can therefore assume that the higher losses in the early autumn were the consequence of the slower reaction of rising groundwater levels surrounding Kozjak Lake in relation to the inflows into the lake system, and consequently of the higher hydraulic gradients of underground run-off from the lake.

The losses from the upper Korana course were much more marked, with the river even running completely dry for several days in the most upstream section of its course. In certain years, such as 2003, the river was dry for a much longer period (57 days), for as many as 36 consecutive days. It was identified that the relation between the discharges measured at the Luketići–Korana hydrological station and the losses was very strong and that the lower the Korana discharge, the higher the share of losses. The average decrease in discharge from the location where the river Korana is formed to the aforementioned Luketići hydrological station was 0.93 m³/s (52% of the losses) in the 1979–1987 period, and 0.94 m³/s (63% of the losses) in the 2015–2020 period, with a total average of 0.93 m³/s (55 %). Looking at the losses in a section of the Korana course further downstream, the Korana–Most profile, the average discharge decrease was 1.12 m³/s (63% of the losses) in the 1979–1987 period and 1.32 m³/s (as much as 89% of the losses) in the 2015–2020 period. The differences in these shares of losses were the result of the different hydrological conditions/inflows of the Korana River during the analyzed periods, though the relation between the discharges and the share of losses did not change during the analyzed periods.

The results obtained from the monitoring of the hydrological conditions in the Plitvice Lakes and the resulting trends of the selected indicators and the changes that they caused showed a strong link between the hydrological conditions and the growth of tufa barriers, as well as the losses of water from the lake system and the Korana River channel. Furthermore, while for the water losses there were no significant differences in the scale of losses during the pre-1990 period and the post-1995 period, huge differences were observed in the barrier growth—the pre-1990 period saw growth of the barrier, but after that period the barriers decreased in size. This indicated the very high sensitivity of the lake system to change, which could in the future be even more marked, due to the expected increasing severity of climate-change impacts. In light of this, in order to protect the lake system, it is necessary to actively monitor its status and any changes that occur. To achieve this, the analyzed example of the Plitvice Lakes can serve as a model for similar locations in the world.