 Definition of the river Gacka springs subcatchment areas on the basis of hydrogeological parameters  Jasmina Lukač Reberski1, Tamara Marković1 and Zoran Nakić2 1 Croatian Geological Survey, Sachsova 2, 10000 Zagreb, Hrvatska; (jlukac@hgi-cgs.hr) 2 Faculty of Mining, Geology and Petroleum Engineering, Pierottieva 6, 10 000 Zagreb, Hrvatska doi: 10.4154/gc.2013.04 Geologia Croatica 66/1 39–53 12 Figs. 4 Tabs. Zagreb 2013 Geologia CroaticaGeologia Croatica Ab sTrA CT The river Gacka springs catchment area is located in the Dinaric karst, which is globally known as the locus typicus, or classical Karst. It is composed of four major and several minor karst springs of different discharge rates. The riv- er Gacka springs are characterised by great discharge and exceptional quality, so the catchment area of the river is indicated in the Water Management Strategy (OFFICIAL GAZETTE NO. 91/08) as an area with strategically im- portant reserves of drinking water for the Republic of Croatia. To determine the hydrogeological characteristics of the subcatchments of this large and complex aquifer system, hydrological and hydrochemical parameters were meas- ured on the main springs. Data collected on the springs were analysed using the recession analysis by the „matching strip“ method, the statistical analysis of a time series of measured data both by autocorrelation (analysis of individ- ual series) and by cross-correlation methods (analysis of interrelationships between time series), multivariate statis- tical analysis (Factor Analysis) of hydrochemical parameters using the software package STATISTICA 6.0 (1998), and geochemical modelling of hydrochemical parameters using the NETPATH computer program. Interpretation of lithological, structural and tectonic characteristics of the rocks, together with tracing data and the applied analytical methods, allowed the springs catchment of the river Gacka to be divided into three subcatchments. The results of this study imply the necessity of a multidisciplinary approach to research. Keywords: hydrogeology, karst spring, hydrology analysis, hydrogeochemistry, river Gacka springs 1. INTrODUCTION The catchment of the river Gacka springs is located in the Croatian part of the Dinaric karst, which represents the karst locus typicus. The catchment consists of four major (Tonko- vić, Majerovo Vrilo, Klanac and Pećina) and several minor springs (Jaz, Marusino Vrilo and Graba) of different dis- charge rates. The Gacka springs are characterised by a high discharge and exceptional water quality, which was why this catchment area was proclaimed by the Water Management Strategy (OFFICIAL GAZETTE NO. 91/08) as one of the strategically important area of drinking water reserves of the Republic of Croatia. The research area is characterised by features of typical karst geomorphology. Research and data collection is difficult due to the heterogeneity of the karst aquifer, as well as the unknown geometry of the voids through which the water flows underground to the springs. In addition, with a size of approximately 500 km2 the river Gacka springs catchment area is categorised as a „big catch- ment” (KENDALL & MCDONNELL, 1998), which further complicates the study. Results of previous hydrological and hydrogeological studies (PAVIČIĆ et al., 2003; LUKAČ RE- BERSKI, 2008) have shown that the river Gacka springs catchment area can be divided into subcatchments of the main springs. Former detailed hydrological analyses were performed on the basis of hydrological data measured only on the profiles of the river Gacka (Bonacci & Andrić, 2008; Lukač Reberski, 2008). Hydrograph analyses performed based on these data are related to the reaction of the entire aquifer system to precipitation, while the characteristics of particular sub-parts of the catchment remain unknown. Therefore, detailed hydrological analyses were performed Geologia Croatica 66/1Geologia Croatica 40 on the most important springs (Tonković, Majerovo Vrilo, Klanac and Pećina), providing new insight on this complex and strategically important karstic hydrogeological system. The investigations and the results presented here are part of an extensive and systematic hydrogeological study described within the authors’ PhD thesis (LUKAČ REBERSKI, 2011). 2. GEOGrAPHIC, CLIMATOLOGIC AND GEOMOrPHOLOGIC CHArACTErIsTICs The investigated area can be described as a hilly karst area with a large karst polje (Gacka polje) surrounded by the Ve­ lebit, Kapela, and Plješivica mountains. The highest altitudes (<1200 m.a.s.l) belong to the central part of the karst hinter- land, whereas the river Gacka springs formed at the lowest altitude, (~450 m.a.s.l). Locally, the relief plays a key role in the definition of climate conditions. High mean values of annual air temperature are the highest in Croatia. There is also a significant difference in the spatial distribution of pre- cipitation. These, are the direct result of the orographic char- acteristics of the study area. The mean annual air tempera- ture ranges between 4 and 9°C, (ZANINOVIĆ et al., 2004) while the annual mean precipitation values range between 1000 and 2500 mm (GAJIĆ­ČAPKA et al., 2003). The re- search area is characterised by typical karst geomorphology i.e. karst features are distinguishable in the relief (karrens, sinkholes, caves, karst poljes, karst springs, estavelles, sinks). 3. GEOLOGICAL, HYDrOGEOLOGICAL AND sTrUCTUrAL-TECTONIC CHArACTErIsTICs The study area is composed predominantly of carbonate rocks with fractured and fracture-cavernous porosity. Over the wider area of the river Gacka springs catchment, various Ju- rassic and Cretaceous limestones and dolomites predomi- nate. The great permeability of the limestone is the direct result of its fragmentation caused by intense tectonics, as well as by the lithological constitution which enables disso- lution processes. Besides the carbonates, Tertiary clastites, i.e. Jelar deposits (BAHUN, 1974) are also abundant in this area, and can be observed overlying the carbonates, some- times as only a few metre thick erosional remains, or up to several hundred metre thick deposits. These deposits are rep- resented by limestone breccias made of unsorted Jurassic, Cretaceous and Palaeogene sediment fragments (VLAHO­ VIĆ et al., 2009). The Jelar formation has specific hydrogeo- logical characteristics and, due to its lithology and spatial distribution, a specific hydrogeological function. The pres- ence of marly-matrix breccias intercalated with marl lenses led to two opposing hydrogeological effects i.e. firstly the unusually well-developed surface and underground karst phe- nomena, and secondly the formation of relatively less per- meable environments (BAHUN & FRITZ, 1975). These de- posits can be observed at the bottom of the polje near the river Gacka spring, as well as in the western part of the catch ment. Fine-grained sediments with intergranular porosity occur in the form of the Quaternary cover in poljes and depressions and they do not have an important hydrogeological influence. Lithology had an essential role in the formation of the hydrogeological properties, especially permeability (LUKAČ RE BERSKI et al., 2009). Hence, based on their permeabil- ity the deposits are divided into: – Very permeable carbonate rocks, which are the most common types, and where limestones of Jurassic and Creta- ceous ages prevail; – Predominantly permeable rocks are represented by other Jurassic and Cretaceous carbonates and Palaeogene clastites, i.e. the Jelar deposits. The Jelar deposits exhibit about 50% less permeability than the Jurassic and Cretaceous carbonates, which surround them (BAHUN & FRITZ, 1975); – Predominantly impermeable rocks – comprising of dolomites of Jurassic and Cretaceous periods; – Impermeable deposits include the Eocene marls oc- curring in a very small region in the eastern part of the in- vestigated area; – Variably permeable deposits are Quaternary deposits with heterogeneous properties, as their permeability varies depending on the thickness and composition. They are rep- resented mostly in the depressions of karst polje. One of the most important, and currently unsolved, karst problems is determining the boundaries and the surface of the catchment, for which classical geological and hydrogeo- logical approaches are necessary but not always sufficient. Hydrological approaches may provide answers to some ques- tions regarding the size of the catchment surface but cannot imply the position of the boundaries (BONACCI, 1995). Un- like a topographic water divide, determination of the ground- water divide in karst terrain depends on many factors. Fur- thermore, the groundwater divide is not constant as its po sition changes depending on the groundwater level (ŽUGAJ, 2000). Therefore a multidisciplinary research approach and data collection is required in order to precisely determine the lo- cation of the divide and the size of the catchment. Several tracings tests were carried out using the Na­fluorescein arti- ficial tracer to determine the apparent velocities and ground- water flow directions as well as the size of the catchment area. Tracing test results, performed between 1957 – 2010., were presented in unpublished technical reports. Recent in- vestigations suggest that the surface of the Gacka catchment up to the Čovići profile extends to 516 km2 (BIONDIĆ et al. 2010). Groundwater­flow directions depend more on the struc- tural relationships than on the lithological characteristics be- cause of the relatively monotonous lithology. The study area is characterized by the NW–SE Dinaric strike of the struc- tures. Tracings showed that the groundwater flow direction in the greater part of the catchment area is parallel to these structures. This is due to the position of the main boundary faults, favourable to the stress direction (20°–45°), and along which the right transcurrent shifts appear, apart from the opening up of spaces and appearance of structures of a pull- apart type (PRELOGOVIĆ, 1989). A fault direction unfa- vourable to the stress direction leads to local compression by closing of the space and possible prevention of water flow (Fig. 1). Lukač Reberski et al.: Definition of the river Gacka springs subcatchment areas on the basis of hydrogeological parameters Geologia Croatica 41 4. HYDrOLOGICAL CHArACTErIsTICs The river Gacka is characterized by favourable hydrological conditions. The flow regime is standardized compared to other karst rivers, and considerable variations of the flow quanti- ties through its river bed do not occur. The ratio of the low- est, medium and highest discharge is 1:4:20, whereas the ratio of the neighbouring Lika river is 1:100:800 (BONACCI, 1987), which shows its vividly torrential character. During the monitoring period (2008 – 2010) the amount of precipitation varied, with 2008/09 being considered as an average year when compared to the long term average, whe- re as 2009/10 was wetter. The average annual discharge dur- ing the monitoring period, (determined from data collected on the Gacka river gauge at Čovići Podgora, Fig. 1), was 13.44 m3s–1 of which 29% is related to the Tonković spring, 26% to the Klanac spring, 24% to the Majerovo Vrilo spring, 15% to the Pećina spring and 6% to the remaining springs (Fig. 2). Recent investigations of the hydrological balance of the catchment indicated the significant retention properties of the karst underground system in the river Gacka catchment Figure 1: Schematic hydrogeological map of the study area, showing sampling sites (position of the study area in detail). 1– very permeable rocks, 2– predominantly permeable rocks, 3–predomi- nantly impermeable rocks, 4– most important faults, 5– extension zone, 6– swallow hole (ponor) and vertical cave (pit), 7– spring: perma- nent and intermittent, 8– groundwater divide, 9– borehole, 10– groundwater connection (tracer test), 11– river, 12– river and rain gauge stations, 13– groundwater velocity determined by tracing tests (cm/s). Geologia Croatica 66/1Geologia Croatica 42 (LUKAČ REBERSKI, 2008). Therefore the water balance of the catchment should be prepared over longer time peri- ods because it provides more reliable discharge coefficient information. Due to the inadequate distribution of rain gauge stations (Fig. 1) in the river Gacka springs catchment, mean an nual precipitation was calculated using the regression equa- tion defined for the research area (GAJIĆ­ČAPKA et al., 2003). This equation also takes into account the precipitation change with altitude as all the rain gauge stations in Lika are below 760 m.a.s.l. and most of the catchment is located at higher altitudes. Mean annual precipitation is 1383 mm while the mean catchment discharge is determined from 30 years of daily discharge data (1972 – 2002), so the mean annual discharge in the profile of the hydrologic gauge station Čovići­Podgora is 14.2 m3/s. The discharge coefficient c can be defined as follows: (1) where Q is the mean discharge (m3s–1), T is duration of the mean discharge (s), P is the amount of precipitation (m), and A is the surface of the catchment (m2). The discharge coefficient of the catchments calculated using the above formula is 0.63, which means that 63% of total precipitation in the catchment infiltrates into the ground and flows towards the river Gacka springs. 5. METHODs AND rEsEArCH TECHNIQUEs TO DETErMINE THE DYNAMIC CHArACTErIsTICs OF THE sPrINGs To determine the dynamic characteristics of individual springs and their catchments, daily water discharge rates were col- lected between 2005 – 2010 from the Meteorological and hydrological service of Croatia (DHMZ). Due to the extre- me ly heterogeneous characteristics of the karst aquifer, spring runoff hydrographs were used because they refer to the total response of the aquifer to precipitation (PADILLA et al., 1994). Analysis of spring hydrographs is both a cheap and very use- ful method in hydrogeological research, especially for karst aquifers. Both the recession parts as well as the complete hydro- graphs were analysed. Recession analysis refers to the anal- ysis of the falling limbs of the runoff hydrographs. The ad- vantage of the recession method is that it does not need prior knowledge of the distribution of the potential and the indi- vidual parameters of the aquifer, although the recession co- efficient directly depends on them. Recession curves mostly represent only part of the dynamic reserves of the aquifer, depending on the initial water level. The total water volume stored inside the aquifer can be much higher. Numerous au- thors have publicised the usefulness of this method (TAL- LAKSEN, 1995; DEWANDEL et al., 2003; BAKALOW- ICZ, 2005). In the case of long recession periods, which a lack of heavy and long precipitation periods, the recession curve could provide good insight into the structural charac- teristics of the aquifer. Both the shape and the characteristics of the recession curve, i.e. the calculated recession coeffi- cient (α) values, imply the hydrogeological characteristics of the aquifer and are dependent on various factors, among which the most important are: porosity type, hydrological conditions, i.e. the water level as well as the rate of water supply from other aquifers. If the α values are larger than 10–2 a rapid drainage system is indicated, with large fractures and channels; when the α value is smaller than 10–2, slow drainage from smaller pores and fractures may occur, i.e. it indicates that the rock matrix has the dominant role in the water flow through the karst underground (KREŠIĆ, 1997). For recession analysis, it is of major importance that the recession period lasts for a long time period, at least for sev- eral months. Such a case has never been registered for the region under study. Precipitation during the recession period results in various shapes of the recession curves for the same spring. Therefore, analysis of a single recession period on a spring is not the most appropriate method for drawing major conclusions regarding the aquifer and groundwater supplies. To overcome this problem the modified „matching strip” meth od was used (POSAVEC et al., 2006) for recession analy­ sis on the monitored springs. This method includes analysis of all recession events of the runoff hydrograph, allowing construction of a main recession curve. In some cases, curves with a lower coefficient of determination describe the reces- sion more appropriately, and this can be easily determined visually (RIGGS, 1968). So the main recession curve, which has the highest coefficient of determination, is not always the best selection. Therefore, the final regression model is selected on the basis of both, experience and the existing software. Recessions with two, three and sometimes even more discharge micro-regimes are quite common in karst aquifers (BONACCI, 1987). Such recession curves can be easily sep- arated. For this purpose software was used for the automatic separation of the main recession curve of the monitored springs (POSAVEC et al., 2010; PARLOV, 2010). The pro- gramme is based on the analysis of the duration curves of mean daily discharge, i.e. on determining the critical dis- charge where the slope of the curve changes abruptly. Analyses of complete hydrographs consist of statistical analysis of time series data (BOX & JENKINS, 1974). This Figure 2: Proportion of discharge from the river Gacka springs between 2008 – 2010. Lukač Reberski et al.: Definition of the river Gacka springs subcatchment areas on the basis of hydrogeological parameters Geologia Croatica 43 method was introduced by MANGIN (1981, 1984) to the karst hydrology and subsequently improved by PADILLA et al. (1994). Numerous authors have used this method for de- fining karst systems (LAROCQUE et al., 1997, 2000; TO- MASKO et al, 2001; PANAGOPOULOS & LAMBRAKIS, 2006; FIORILLO & DOGLIONI, 2010; STROJ, 2010; TER- ZIĆ et al., 2011)., Autocorrelation and cross­correlation meth ods are used for the time series analysis of the springs, which is the aim of the current research (DAVIS, 2002; BOR- RADAILE, 2003). The characteristics of the springs catch- ment, the relationships between springs and the dependence between the precipitation in the catchment and the spring discharge were all determined using these methods. The autocorrelation function represents the linear rela- tionship of the successive data values inside a time series, dependant on their time distance (BOX & JENKINS, 1974; TERZIĆ et al., 2011). The autocorrelation method is based on the intercorrelation of data from a time series, with a suc- cessive increase of time distances, whereby the degree of similarity between the data for each step is calculated. Au- tocorrelation for the distance τ corresponds to the covariance of all measurements xt and the measurements with a time distance xt + τ according to the equation: (2) where cov τ is the covariance for a time distance, x is the time series, n is the number of measurements in a time se- ries, τ is the time distance between two measurements, is the average value of the sample. The autocorrelation coefficient (rx) ranges from a max- imum of 1 to the minimum of –1. The autocorrelation value of 1 means that the compared time series are identical. The time needed for the decrease of the correlation coefficient below 0.2 is called the „memory effect” MANGIN (1984), and reflects the duration of the reaction of the system to the input signal. Interpretation of the autocorrelogram should be undertaken with care, because the shape of the autocorrelo- gram depends, not only on the characteristics of the karst system, but also on the duration, frequency and the intensity of the rain (GRASSO & JEANINE, 1994; EISENLOHR et al., 1997). To establish the places of pronounced similarities i.e. similarities between individual data, it is possible to compare two different time series in the same way one time series is compared with itself for a particular time lag. From such comparisons it is possible to obtain two kinds of informa- tion: the strength of the connection between these two time series and the time delay between them at the places of max- imum similarity. The cross-correlation equation is similar to the autocorrelation coefficient equation. If two time series are marked as variables X and Y, and n is the number of pairs which are compared in one step (k) of the cross-correlation, than the cross­correlation coefficient is as follows: (3) 6. METHODs UsED FOr DEFINING THE HYDrOGEOCHEMICAL CHArACTErIsTICs OF sPrING WATErs Between February 2008 and November 2010, water samples were collected from four springs (Tonković, Pećina, Maje- rovo Vrilo and Klanac). Prior to sampling, the following pa- rameters of the sample waters were measured ‘’in situ’’ by probes from the WTW company: electrolytic conductivity (EC), temperature (T), pH and oxygen content. Water samples (500ml) were collected for measurements of basic cations, anions and alkalinity. Samples (100ml), were also collected for measurements of stable isotopes (δ18O and δD). The sam- ples were collected under different hydrological con ditions, with a total of 16 samples per spring for basic cations and 25 samples for stable isotope measurements. Basic cations and anions were measured in the Hydrogeochemical labora- tory, Department of Hydrogeology and Engineering Geol- ogy of the Croatian Geological Survey. Con centrations of sodium, potassium, magnesium and calcium were measured using the atomic absorber type AAS Analyst 700 by Perkin Elmer. Nitrate, orthophosphate, sulphate and chloride ions were measured using the ion chromatograph by LabAlliance. Bicarbonate ions were measured by titration. From the che- mi cal composition of the water, interactions between the rock-mass and groundwater can be determined, as can the underground retention time of the water. Besides the natural processes, the water composition is also influenced by pol- lution, especially by nitrates found in artificial and natural fertilizers. Understanding of the reactions between the water and environmental substances introduced to it allows antic- ipation of the dangers/harmful influences of pollutants on water quality, and hence the influence on the human health and/or the health of the ecosystem. The ratios of the stable isotopes deuterium (dD) and ox- ygen (d18O) in the sampled waters were measured using the mass spectrometar by Finnigan DELTA plus XP i DELTA plus in the JOANNEUM RESEARCH RESOURCES, Institute of Water, Energy and Sustainability, Department of Isotope Hydrology and Environmental Analytics in Graz (Austria) and in the Physics department of the School of Medicine in Rijeka (Croatia). These ratios are used in hydrogeological research to determine the recharge area and the hydrodyna- mic conditions, i.e. changes in water velocity in discrete parts of the aquifer and the hydraulic connection between separate aquifers (DOCTOR et al., 2000; VANDENSCHRICKA et al., 2002). In nature the differences in the ratio can be seen in the phase transitions in the water media. Therefore, vari- ations of the isotope ratios in the groundwater are the result of seasonal changes, the altitude of the recharge area, prox- imity to the sea, intensive rain showers and isotope exchange with the rock-mass. The linear dependence between the sta- ble isotopes oxygen-18 and deuterium in precipitation was recognized by FRIEDMAN (1953), while CRAIG (1961) defined this dependence with the following formula: dD = 8 · d18O + 10 (4) Geologia Croatica 66/1Geologia Croatica 44 which he named global meteoric water line (GMWL). The GMWL is calculated for the mean ratios of oxygen and hy- drogen isotopes from precipitation around the world. How- ever, the oxygen and hydrogen isotope ratios of some regions can differ from the global ones due to different meteorolog- ical conditions in the region. In this case, the local meteoric water line (LMWL) can be calculated. The physico-chemical indicators which were continu- ously collected during fieldwork, as well as hydrochemical measurements, were geochemically modelled using the NET PATH computer programme and statistical multivariate Factor Analysis using the STATISTICA software package. NETPATH was developed by PLUMMER et al. (1994) from the USGS and is based on the geochemical model of mass balance. It consists of a series of smaller programmes which, when put together, enable the input and management of che- mical and isotope data from the analysed water. Comparison of the determinants between samples was undertaken by appli- ca tion of the R-mode Factor Analysis, using the STATISTI CA 6.0 programme package (1998), and showed geochemical pro- cesses were the dominating ones determining the che mical composition of the water,. Factor analysis was performed on the following indicators: T, pH, Ca2+, Mg2+, nMg/nCa, Na+, K+, HCO3 –, Cl–, SO4 2–, NO3 –, PO4 3–-P, SIcal, SIdol, logpCO2. 7. rEsULTs AND DIsCUssION 7.1. recession analyses The recession coefficients of the monitored springs were ob- tained using the „matching strip” method for the analysis of the falling limbs of the runoff hydrographs (Fig. 3). Since only the exponential regression model enables evaluation of the recession coefficient on which comparison of the prop- erties of the aquifer are based,, this model was selected as being the most appropriate for comparison of the characteris- tics of the catchments of the monitored springs. For the dis- charge time series measured at Majerovo Vrilo, the most ap- propriate solution for the recession analysis was dividing the main recession curve into two parts at the point where the critical discharge rate was 6,4 m3/s. The value of the reces- sion coefficient of the first part of the curve is α1 = 0,247 day–1, and the second α2 = 0,029 day–1. For other springs, division using three discharge micro-regimes was selected. At the Klanac spring the discharge regime changes at the value of 4,13 m3/s. This implies that high waters are best described by the recession coefficient of α1 = 0,044. Changes in the runoff regime also appear at the critical discharge value of 1.14 m3/s. Hence the recession coefficient of α2 = 0,032 was calculated for the medium waters. The recession coefficient of low waters is α2 = 0,046. In the case of the Pećina spring, the main recession curve was divided at the critical discharge value of 6.00 m3/s, and the high waters are best described by the recession coefficient of α =0,115. The medium waters recession coefficient value is α2 = 0,047 lasting until the crit- ical discharge value of 0,65 m3/s. The value of the recession coefficient for low waters is α3 = 0,07. Discharge of the aqui- fer at the Pećina spring during the middle water regime is not uniform (Fig. 3), suggesting that the regression model selected on the coefficient of determination is not significant. At the Tonković spring, changes in the discharge regime appear at values of 4.20 m3/s and 2,70 m3/s. The recession coefficient defined by the automatic recession analysis is α1 = 0,023 for high waters, α2 = 0,014 for medium waters and α3 = 0,009 for low waters. As at the Pećina spring, dur- ing the high water regime, significant differences can be de- termined in the drainage speed during individual recession periods at the Tonković spring. This was the reason why it was not possible to analyse the high water recession curve with more certainty. Figure 3: Recession analyses at the observed springs. Lukač Reberski et al.: Definition of the river Gacka springs subcatchment areas on the basis of hydrogeological parameters Geologia Croatica 45 Recession analysis (Fig. 3, Table 1) clearly showed that the Pećina and Majerovo Vrilo springs have significant dif- ferences between the recession coefficients of the quick and the base flow, when compared to the Klanac and Tonković springs. The highest runoff differences between the base and quick flow occurs at the Majerovo Vrilo spring with the dif- ference in the recession coefficient between the base and the quick flow being of one order of magnitude. Furthermore, as seen on the recession curve graphs of the springs, the Pe- ći na and Tonković springs have significant differences in the types of drainage of the fissured systems during different quick flow recession events. At the Pećina spring, during in- dividual recession periods, the runoff decreased from 5 m3/s to 1 m3/s in only 4 – 5 days, while during other recession periods for the same runoff more days are needed. At the Tonković spring, significant differences in the discharge rates during different recession periods of the high water regimes can also be determined. While in some cases of recession, run- off decreases from 7 to 4 m3/s in only 2 days, sometimes 20 days were needed for the same decrease. This variation in the runoff properties at the same spring could be the result of: the varying intensity and spatial distribution of precipitation in the catchment, the land cover, differences in evapo trans- piration rates, and the heterogeneity of the aquifer, which di- rectly influences the dynamics of the runoff. At the Klanac and Majerovo Vrilo springs, no such significant differences in the runoff during different recession events were determined. At the Tonković spring the autocorrelation function (AKF) has a steeper slope for the first 20 days. The memory effect, i.e. the average duration of the springs reaction to pre- cipitation is 62 days. Similarly, the slope of the function is steeper during the first 22 days at the Klanac spring and the memory effect is 63 days. At the Majerovo Vrilo and Pećina springs, the shape of the autocorrelation functions are some- what different. At the Majerovo Vrilo spring the steepest slo pe of the function is in the first 4 days, while the response of the catchment, from which the spring is recharged, is 26 days. At the Pećina spring the fast discharge lasts 9 days on average, and the memory effect of the spring is 33 days. Steep slopes of the autocorrelation function (AKF) indicate rapid infiltration of precipitation and faster drainage of the underground system through the well developed network of karst channels. Its more gentle slopes reflect slower drainage of the fissured porosity aquifer area where discharge is sim- ilar to the intergranular aquifers. It was established, after the statistical analysis of the dis- charge time series using the autocorrelation method, that all four investigated springs drain water from aquifers with sim- ilar hydrogeological properties. This is evident from the sim- ilarities of their autocorrelation functions (AFK) (Fig. 4). Nevertheless, at the Pećina and Majerovo Vrilo springs the AKF have a steeper slope at the beginning and the values of the memory effect are lower than at the Klanac and Tonković springs. This is the result of the better developed underground karst channel network of the Pećina and Majerovo Vrilo aq- uifers, which is in accordance with the results of the reces- sion analysis. When interpreting AKF, it should be known that they relate to the average properties of the analyzed time series. From the cross-correlation functions (KKF) of discharge time series of the different springs, the connection and the interrelation of the springs can be recognized (Fig. 5). Sym- metrical KKF values of all springs indicate that the springs are independent of each other and depend on some other va- riable. In this case it refers mostly to precipitation in the catch ment of the monitored springs, because their quantity and time and space distribution, together with the hydrogeo- logical aquifer properties, influence the runoff dynamics at the springs. It is supposed, based on the high values of the maximum cross­correlation coefficients (max rxy), which range from 0,77 – 0,91, that all springs drain aquifers from catchments of similar hydrogeological characteristics and Table 1: The values of the recession coefficients at observed springs. Recession coefficient Pećina Majerovo Vrelo Tonković Klanac α1 0.115 0.275 0.023 0.044 α2 0.047 0.029 0.014 0.032 α3 0.069 0.009 0.046 Table 2: The values of basic statistical parameters for the flow rate. Qmin (m3/s) Qmax (m3/s) Qx (m3/s) Qσ (m3/s) cv (%) Klanac 0.017 16.70 3.466 3.253 93.85 Pećina 0.002 6.84 1.925 1.694 88.01 Tonković 1.680 11.40 3.943 1.980 50.22 Majerovo Vrelo 1.000 21.03 3.264 2.263 69.32 In cases of extreme differences in the recession charac- teristics of the aquifer, it is not possible to determine the re- cession curve with a high value of the coefficient of deter- mination R2. Nevertheless, the „matching strip” method is still considered to be the more appropriate solution for the determination of the recession properties of the aquifer com- pared to analysing individual recession events. For the anal- ysis of individual recession events, it is not possible to de- termine how much the runoff regimes of different recession events at the same spring differ from each other. 7.2 Analysis of the complete hydrographs From the basic statistical water flow parameters presented in Table 2 it can be seen that the mean discharge values (Q) of the Klanac, Tonković and Majerovo Vrilo springs are very si­ milar, while at the Pećina spring it is almost twice as low. The values of the variation coefficient (cv), which represents the ratio between the standard deviation (σ) and the mean value ( ) of the sample, indicate that the highest discharge varia- tions occur at the Klanac, and lowest at the Tonković springs. Geologia Croatica 66/1Geologia Croatica 46 with a similar precipitation volume. The delayed response to rainfall of the Klanac and Tonković springs if compared to the Majerovo Vrilo spring, ranges from 24 – 36 hours. The delay at the Pećina spring is between 12 – 24 hours. The Ma- jerovo Vrilo spring reacts up to 12 hours prior to the Pećina spring, the Tonković spring up to 12 hours prior to the Kla- nac spring. The highest cross­correlation coefficients are those of the Klanac and Tonković springs, with the lowest at Pećina and Majerovo Vrilo, which is expected since these two springs have the most distant catchments. The maximum discharge and precipitation cross-corre- lation coefficients max rxy (Fig. 6) are not high (0.3 – 0.46). The reason for the low values is that the discharge dynamics at the springs does not solely depend on the quantity of pre- cipitation in the catchment. It also depends on the season, i.e. air temperature during precipitation, the vegetation and the hydrogeological properties, including the degree of karsti- fication, the existence of cover deposits, their thickness and other factors. Nevertheless, despite the low cross-correlation coefficients the connection between discharge and precipita- tion is clearly pronounced. At the Pećina spring, the maxi- mum cross­correlation coefficient was registered during the three day delay. This means that the average reaction delay of this spring is three days. At the Majerovo Vrilo spring, the reaction to precipitation is the same as for the Pećina spring, but the max rxy is somewhat higher than the other springs. The cross­correlation coefficients of the Klanac and Tonković springs have two maxima which means that different pre- cipitation events cause variation in delay. Hence the Klanac spring reacts in 4 to 10 days, and Tonković spring 5 to 11 days after the onset of precipitation. Various delay times at these springs could be the result of variation in infiltration of precipitation in the karst underground over both time and place in the catchment. Rapid reaction of a spring to precip- itation indicates fast infiltration in the epikarst zone and fur- ther into the underground channel system, which enables ra- pid transport. Another reason could be a smaller catchment surface which makes the path the infiltrated water must travel from the surface to the spring significantly shorter. After a time period of 40 days the function gradually starts to grow again as seen on the cross-correlograms. This is in accord- ance with the autocorrelation functions of the springs, i.e. it is the result of repeated incoming water waves, especially in the wetter part of the year. All hydrograph statistical analysis, as well as recession hydrograph analysis used in this study, indicate faster reac- tions to rainfall by the Majerovo Vrilo and Pećina compared to the Tonković and Klanac springs. This implies a better developed underground channel network of the aquifer sys- tems of the Majerovo Vrilo and Pećina springs. The Klanac Figure 4: Autocorellograms of spring daily discharges. Figure 5: Cross-corellograms of discharge series. Lukač Reberski et al.: Definition of the river Gacka springs subcatchment areas on the basis of hydrogeological parameters Geologia Croatica 47 and Tonković springs aquifer systems show better retention properties and lower permeability due to the network of smal ler fractures within the aquifer. 7.3 The isotopic composition of the waters The LMWL was calculated for the research area by MANDIĆ et al. (2008) (Fig. 7) who emphasised that there is a need for several more years of measurements in order to determine a more precise value. The LMWL equation is as follows: dD = 6.8·d18O + 1.5 (5) The measured stable isotope oxygen-18 and deuterium values of the spring waters are scattered around the LMWL and GMWL values, clearly indicating that in the study area the groundwater is recharged from precipitation (Fig. 7). It was shown by MANDIĆ et al. (2008) that the average deu- terium excess in the investigated area is about 12 ‰. Such values are evidence of the influence of the Mediterranean type of climate, especially during the winter months, but the main influence is from the continental type of climate. No significant fluctuations in oxygen­18 values (Fig. 8) were observed in the spring waters of the river Gacka catch- ment which is in agreement with the observations of MAN- DIĆ et al. (2008). The lowest negative d18O values were mea sured at Majerovo Vrilo, followed by those at the Tonko­ vić spring. The Pećina spring has the highest positive values. These values imply the highest altitude of the recharge area occurs at the Majerovo Vrilo and the lowest at the Pećina spring. For Europe and Africa, the δ18O isotope gradient ran- ges from –0.15 to –0.5 ‰ for every 100 m increase in alti- tude (CLARK & FRITZ, 1997). In order to calculate the al- titude of the recharge area, the isotope gradient of –0,25 ‰ was selected and the δ18O values measured in spring and autumn, as suggested by CLARK & FRITZ (1997). The al- titude differences obtained range from 140 to 260 m. CLARK & FRITZ, (1997), determined that the depletion for 18O var- ies between about –0.15 to –0.5 ‰ per 100 m rise in altitude. The calculated altitude values of the dominant recharge area of the springs are as follows: • Majerovo Vrilo – 1020 m a.s.l., • Tonković – 760 m a.s.l. • Pećina – 620 m a.s.l. 7.4. The quality of the spring waters The results of the measurements of the water quality indica- tors are given in Table 3. Based on their basic chemical composition the spring waters of Majerovo Vrilo belong to the Ca­HCO3 to CaMg- HCO3-type of water, i.e. to the hydrochemical facies (Fig. 9). Spring waters of Tonković, Klanac and Pećina belong to the Ca-HCO3 type of water. These water types are the result of the dissolution of carbonate minerals (calcite and dolo- mite) which are the main constituents of the catchment area of the springs. Figure 6: Crosscorellograms of springs which represent the relationship between rainfall and discharge. Figure 7: The ratio of oxygen-18 and deuterium stable isotopes in spring waters of the Gacka catchment. Geologia Croatica 66/1Geologia Croatica 48 The temperature of the spring waters, was measured pe- riodically during the research period, and varies from 8.4 to 10.0 °C. This is generally in accordance with the mean an- nual air temperature of the recharge area of the spring. The lowest average measured temperatures were recorded at Ma- jerovo Vrilo, which is in accordance with the data suggest- ing that this spring is recharged from the highest altitudes, where the temperatures are the lowest. Temperatures meas- ured at the Tonković and Pećina springs are also compatible with the calculated altitudes of their recharge area. The oxygen content measured in the spring waters, as well as the model determined carbonate saturation level, is advantageous, from the water quality protection standing point. A high oxygen concentration was measured in the wa- ters which enables oxidation processes useful as water au- topurification reactions. Spring waters are mostly saturated with calcite, which enables deposition of any existing metals in the water due to their affinity with carbonates. Hence, in saturated conditions these reactions could be initiated, which also has a water autopurification role. The measured nitrate concentrations are below the Croa- tian drinking water standard. (OFFICIAL GAZETTE NO. 47/08) and range from 1.3 to 5.2 mg/l with small seasonal oscillations. The highest concentration was measured at the Pećina spring (Table 3) in the autumn and winter. This is the consequence of both anthropogenic influences and degrada- tion of the plant material in the highly karstified aquifer of this spring. Ammonium and orthophosphate concentrations are very low in all the springs and sometimes they are below the de- tection limit. Figure 8: The distribution of stable isotope oxygen-18 in the spring waters. Table 3: Results of measurements of hydrochemical indicators. Sp rin g T EC pH Ca2+ Mg2+ Na+ K+ HCO3 – Cl– SO4 2– NO3 – PO4 3–-P O2 (°C) (μS/cm) (mg/l) (mg/l) (mg/l) (mg/l) (mg/l) (mg/l) (mg/l) (mg/l) (mg/l) (mg/l) M aj er ov o Vr el o Min 8.8 307 7.22 54.8 2.1 0.8 0.4 204 1.3 3.4 2.1 0.02 8.2 Max 9.6 461 7.69 88.4 26.5 4.5 0.9 312 4.0 5.8 4.3 0.04 11.7 Avg 9.1 437 7.40 71.9 8.7 2.3 0.7 258 2.5 4.5 3.5 0.03 9.9 Kl an ac Min 9.0 440 7.11 62.9 4.4 1.0 0.4 235 1.7 4.3 1.6 0.01 7.2 Max 10.0 493 7.63 138.0 13.8 4.6 1.1 428 4.9 11.4 3.7 0.05 14.4 Avg 9.5 474 7.36 86.5 7.0 3.0 0.7 284 3.4 6.4 3.0 0.03 9.3 To nk ov ić Min 9.5 473 7.10 62.8 4.3 2.5 0.4 235 2.6 4.1 1.3 0.01 6.1 Max 9.9 532 7.60 124.0 14.6 6.9 1.0 398 14.1 12.9 3.9 0.04 11.7 Avg 9.7 495 7.35 86.6 7.2 4.1 0.7 287 5.5 7.3 3.1 0.03 9.0 Pe ći na Min 8.4 442 6.95 60.4 3.0 2.2 0.5 218 2.9 2.9 2.2 0.01 7.2 Max 9.9 503 7.96 138.0 13.4 6.2 1.0 426 5.0 10.9 5.2 0.05 11.2 Avg 9.3 468 7.55 86.0 6.1 4.1 0.7 283 4.2 5.2 3.9 0.03 9.2 Lukač Reberski et al.: Definition of the river Gacka springs subcatchment areas on the basis of hydrogeological parameters Geologia Croatica 49 Sulphate concentrations range from 3.4 – 12.9 mg/l and are below levels for the Croatian drinking water standard (OFFICIAL GAZETTE NO. 47/08). Concentrations of chlo- ride range from 1.3 – 14.1 mg/l with the highest concentra- tion being observed at the Tonković spring in March (14.1 mg/l) and April (8.9 mg/l). Furthermore, concentrations of chloride at other springs are higher in March and April, than during other months. These high concentrations are most like ly to be a consequence of de-icing of the roads, because the springs are located near them. Waters of the river Gacka springs are still characterised by their excellent quality and are classified as strategically important water reserves of the Republic of Croatia (OFFI- CIAL GAZETTE NO. 91/08). Hence, it is important to mon- itor and to promptly react to existing and potential pollution sources in the catchment. Taking into account previous hy- drogeochemical studies carried out in this area (LUKAČ RE BERSKI, 2008; LUKAČ REBERSKI et al., 2009), it is evident that, despite the good protection and quality, the an- alysed water quality indicators show an upward trend due to the anthropogenic influence in the spring catchments. 7.5. statistical analysis of the geochemical indicators of spring waters Five factors were selected based on the results of the Factor Analysis i.e. the factor loadings (Table 4). The first is the lithogenic factor and refers to the influence of the partial pres sure on the calcite and dolomite saturation indices, i.e. the influence on the dissolution/precipitation of limestone and dolomite. The second factor is also a lithogenic factor, and represents the influence of dolomite dissolution on the water chemistry. The third factor is anthropogenic; it shows the influence of the wash out of the terrain and epikarst zone on the chemistry of the water. The fourth factor is lithogenic and represents the dissolution of limestone. The fifth factor is anthropogenic and indicates the influence of salting of roads on the chemical composition of the water. 7.6. The calculations of mixing ratios Mixing models for individual springs were provided using the NETPATH programme. The water ratios from the Maje- rovo Vrilo and Tonković sprin catchments in the Klanac spring were determined in order to establish if this spring shares the catchment area with the others. The water ratios vary depend- ing on the hydrological conditions. During high and medium waters the Klanac spring mostly shares the catchment area with the Tonković spring (Fig. 10). During low water peri- ods it shares its catchment area equally with Majerovo Vrilo and the Tonković springs. It should be noted that these val- ues are valid for the time of sampling. For a more detailed analysis of their interrelationships samples should be col- lected on a daily basis throughout the hydrological year. 7.7. The subdivision of the catchment The catchment of the river Gacka springs was divided into three subcatchments based on the results of the investiga- tions, the interpretation of the lithological, structural and hy- drogeological characteristics of the aquifer system and the hydrogeochemical properties of the spring waters, (Fig. 11). The Tonković and Klanac springs, which are about 50 m away from each other, drain the same catchment as shown by the research results. Hence, it was not possible to estab- lish the groundwater divide between these two springs. The Figure 9: Piper diagram showing the hydrochemical facies of the observed springs (e.g. KLAN-11-08 means SPRING-month-year). Table 4: Factor loadings. Factor score > 0,5 (varimax nom.) Factor 1 Fakcor 2 Factor 3 Factor 4 Factor 5 T (°C) 0.00 –0.12 –0.63 –0.16 0.35 pH 0.17 0.09 0.12 –0.36 0.09 Ca2+ (mg/l) 0.33 0.33 –0.02 0.83 0.22 Mg2+ (mg/l) 0.07 –0.98 0.00 0.05 0.03 nMg/nCa –0,04 –0.97 0.02 –0.17 –0.06 Na+ (mg/l) 0.07 0.08 0.80 0.15 –0.01 K+ (mg/l) 0.25 0.05 0.13 –0.06 0.70 HCO3 – (mg/l) 0.35 0.02 0.07 0.87 0.19 Cl– (mg/l) –0.04 –0.01 –0.21 0.15 0.81 SO4 2– (mg/l) –0.09 0.02 –0.80 –0.02 –0.37 NO3 – (mg/l) 0.43 0.08 –0.18 –0.05 –0.16 PO4 3–– P (mg/l) –0.08 0.00 0.13 0.21 0.38 SI cal 0.90 0.06 0.16 0.25 0.17 SI dol 0.90 –0.16 0.07 0.03 0.17 Log pCO2 –0.53 0.18 0.16 0.63 0.31 Expl. Var 2.42 2.12 1.86 2.18 1.82 Prp. Totl 0.16 0.14 0.12 0.15 0.12 Geologia Croatica 66/1Geologia Croatica 50 groundwater divide between the catchment of the Pećina spring and the Tonković and Klanac springs is not a line, but represents a recharge zone which is determined between the se springs and is common to both catchments. Most of the Tonković and Pećina springs catchment area drains the Tonković spring. During high water levels a connection ex- ists between the sink in Kozjan with the Pećina spring, as evidenced by tracing test results. Furthermore, these test re- sults indicate that the Pećina spring has a larger catchment surface than expected based on the mean annual discharge values. Considering the complex structural-tectonic relation- ships, which are a reflection of many phases of karstification and tectonism, particularly the neo tectonic dynamics (PRE- LOGOVIĆ, 1989), it is difficult to determine the cause of this discrepan with a high degree of certainty. There are several possible explanations. One could be that during pe- riods of heavy precipitation, the amount of rainfall exceeds the infiltration rate and the water level rises in the epikarst Figure 11: River Gacka springs subcatchment ar- eas. 1-subcatchment areas: A Majerovo Vrilo, B Tonković and Klanac springs, C Pećina spring, 2-common catchment area of adjacent subcatch- ments, 3-spring, 4-assumed groundwater flow di- rection, 5-Gacka springs catchment groundwater divide, 6-approximate divide of common areas. Figure 10: Plot of the calculated mixing ratios of waters from the Tonković and Majerovo Vrilo springs in the water from Klanac spring. Lukač Reberski et al.: Definition of the river Gacka springs subcatchment areas on the basis of hydrogeological parameters Geologia Croatica 51 zone. This causes water overflow into the preferred flows near the surface, and possibly establishing a connection with the Pećina spring, probably due to the infiltration through the shallower parts of the aquifer. Another reason could be water flowing over a supposed hydrogeological barrier in the spring area during high waters. Furthermore, underground discharge at the Pećina spring is not excluded due to the frac- tured nature of the area surrounding the spring. This could be the case especially during low waters and could occur without any obvious surface runoff. In that case a larger catch ment surface is drained without any visible evidence on the surface. This could also be in accordance with the re- cently postulated theory (BONACCI & ANDRIĆ, 2008) which assumed that part of the waters of the river Lika dis- appears and flows underground to the river Gacka catchment, draining downstream in the Pećina spring. Three permanent springs are located downstream of the Tonković spring: Jaz, Marusino Vrilo and Graba. They are characterised by small fluctuations in discharge under different hydrological condi- tions which indicate the existence of a barrier in their hin- terland which slows down the groundwater and prevents sig- nificant water level fluctuations. A common zone was also defined between the Majerovo Vrilo catchment and the Klanac and Tonković springs. Water tracings performed in Vrhovinsko polje and at Trnavac con- firmed an underground connection common to both catch ments. Furthermore, the hydrogeochemical mixing model showed that during high water levels the groundwater flows to the Klanac spring from the same part of the catchment as to the Tonković spring. Under low water conditions, at the Klanac spring there is an increased ratio of groundwater from the Majerovo Vrilo subcatchment. This is probably due to hydraulic conditions which prevent water recharge from the Majerovo Vrilo catch- ment during the high water regime. When the water level drops, water recharge from Majerovo Vrilo (Fig. 12) is enabled. It is supposed that the water from the Klanac catchment does not discharge towards the Majerovo Vrilo. 8. CONCLUsION Applied research techniques have enabled the recognition of subcatchments in, parts of the large karst aquifer system of the river Gacka springs and the determination of their dy- namic characteristics, as well as their interference. All the former hydrological analyses were performed based on hy- drological da measured only on the river Gacka profiles. Such data were referred to the entire system, while the char- acteristics of certain parts of the catchment still remain un- known. Therefore, in the cases of such large catchment areas, it is very useful to divide it into smaller units or subcatch- ments. Results presented here show that there are significant differences in hydrogeological characteristics of individual parts of the large river Gacka springs catchment area. Based on interpretation of the lithological and structural-tectonic characteristics of the area, the tracing test data, and the ap- plied research methods, the river Gacka catchment is divided into three subcatchments: the Majerovo Vrilo, the Tonković and Klanac, and the Pećine subcatchments. The boundaries between them are not presented as lines, but were determin ed as zones which are common to the neighbouring subcatch- ments (Fig. 11). Statistical analysis, analyses of the recession parts of the hydrograph, together with hydrogeochemical analyses, all indicate a better development of the groundwater flow paths of the aquifer systems in the Majerovo Vrilo and Pećina springs catchments when compared to the Tonković and Kla- nac springs. This is the reason why the Majerovo Vrilo and Pećina spring catchments are more sensitive to pollution. Due to their well developed network of underground water flow and the short subsurface residence time of the water, they have a low water autopurification capacity. The Tonković spring during the quick flows, showed significant variation in drainage during different recession events, as observed from recession analysis. This led to the conclusion that parts of this subcatchment are even more sensitive to pollution. The autocorrelati and cross-correlation time series anal- ysis and the matching strip method are very useful methods for the hydrogeological research of karst aquifers. They en- able collection of data regarding the average aquifer charac- teristics which can be used f interpretation of the hydroge- ological properties of the catchment, such as the structural properties, type of porosity, development of the channel net- work, dynamics of the aquifer system, sensitivity to pollu- Figure 12: A conceptual model of the river Gacka springs catchment area. Geologia Croatica 66/1Geologia Croatica 52 tion and others. Analysis of the recession hydrograph using the matching strip method proved to be very useful because, in addition to calculation of the recession coefficient, the diagram can be used to establish the time period of various recession events, since all recession events of the time series are presented on the diagram. Hence, such recession analy- sis is more reliable than cases where several recession events are considered individually. The water of the river Gacka springs is still of excellent quality. Hydrogeochemical investigations imply a trend of growing concentrations of anthropogenic indicators (nitrates, orthophosphates, sulphates and chlorides) in the spring’s catch ments. Furthermore, the dominant processes which in- fluence the composition of the spring waters were determin ed, together with the average recharge altitudes of the springs, and using a geochemical mixing model, the interference be- tween the catchments was also established. Results presented here point out the necessity for a mul- tidisciplinary research approach, especially in order to de- termine the hydrogeological characteristics of complex aq- uifer systems like the river Gacka springs catchment. Due to the high importance of this catchment, which was proclaimed to be an area of strategic drinking water reserves for the Re- public of Croatia, the collected data are useful for all future development plans, both for groundwater exploitation as well as for protection from possible contamination. rEFErENCEs BAHUN, S. (1974): Tektogeneza Velebita i postanak Jelar naslaga [The tectogenesis of the Mt. Velebit and the formation of Jelar deposits – in Croatian].– Geol. vjesnik, 27, 35–51. BAHUN, S. & FRITZ, F. (1975): Hidrogeološke specifičnosti Jelarna- slaga [Hydrogeological specificities of the Jelar breccias – in Croa- tian].– Geol. vjesnik, 28, 345–355. BAKALOWICZ, M. (2005): Karst groundwater: a challenge for new resources.– Hydrogeology Journal 13, 148–160. doi:10.1007/ s10040 -004-0402-9. BIONDIĆ, R., BIONDIĆ, B. & MEAŠKI, H. (2010): Novelacija gran- ica zona sanitarne zaštite izvorišta Gacke [Gacka spring sanitary protection zone novelation – in Croatian]. Geotehnical faculty, Uni- versity of Zagreb, 115 p. BONACCI, O. (1987): Karst Hydrology with Special Reference to the Dinaric Karst.– Springer­Verlag­Berlin­Heidelberg, 9/4, 328–338. BONACCI, O. (1995): Međuzavisnost geologije, hidrogeologije i hidro- logije pri rješavanju hidrotehničkih problema. [The interdependence of geology, hydrogeology and hydrology in solving hydraulic proble – in Croatian].– First Croatian geological congress, Book of Pro- ceedings, 1, 105–108. BONACCI, O. & ANDRIĆ, I. (2008): Sinking karst rivers hydrology: case of the Lika and Gacka/Croatia.– Acta Carsologica, 37/2, 185– 196. BOX, G.E.P. & JENKINS, G.M. (1974): Time Series Analysis: Forecast- ing and Control.–Holden Day, San Francisco, 575 p. CLARK, I. & FRITZ, P. (1997): Environmental Isotopes in Hydrogeol- ogy.– CRC Press, New York, 20, 174–188. CRAIG, H. (1961): Isotope variations in meteoric waters.– Science, 133, 1702–1703. DAVIS, J.C. (2002): Statistics and Data Analysis in Geology. Cansas geo lo gical survey.– The University of Cansas, third edition, USA, 638 p. DEWANDEL, B., LACHASSAGNEB, P., BAKALOWICZ, M., WENGB, P. & AL-MALKI, A. (2003): Evaluation of aquifer thickness by ana lys ing recession hydrographs. Application to the Evaluation Oman ophiolite hard­rock aquifer.– J. Hydrol., 274, 248–269. doi: 10.1016/S0022-1694(02)00418-3. DOCTOR, D. H., LOJEN, S. & HORVAT, M. (2000): A stable isotope investigation of the classical Karst aquifer: evaluating karst ground- water components for water quality preservation.– Acta carsologi- ca, 29/1, 5, 79–92. EISENLOHR, L., BOUZELBOUDJEN, M., KIRALY, L. & ROSSIER, Y. (1997): Numerical versus statistical modelling of natural respon­ se of a karst hydrogeological system.– J. Hydrol. 202, 244–262. FIORILLO, F. & DOGLIONI, A. (2010): The relation between karst spring discharge and rainfall by cross-correlation analysis (Campania, south- ern Italy).– Hydrogeology Journal, 18/8, 1881–1895. doi:10.1007/ s10040-010-0666-1. FRIEDMAN, I. (1953): Deuterium content of natural waters and other substances.– Geochim. et. Cosmochim. Acta, 4, 89–103. GAJIĆ­ČAPKA, M., PERĆEC TADIĆ, M. & PATARČIĆ, M. (2003): Digitalna godišnja oborinska karta Hrvatske. [A Digital Annual Pre- cipitation Map of Croatia – in Croatian] Hrvatski meteorološki ča­ so pis (Croatian Meteorological Journal), 38, 21–33. GRASSO, D.A. & JEANNIN, P-Y. (1994): Etude critique des méthodes d’analyse de la réponse globale des systèmes karstiques. Applica- tion au site de Bure (JU, Suisse).– Bulletin d’Hydrogéologie de l’Uni versité de Neuchâtel, 13, 87–113. KENDALL, C. & MCDONNELL, J. J. (1998): Isotope Tracers in Catch- ment Hydrology.–Elsevier Science B.V., Amsterdam, 839 p. KREŠIĆ, N. (1997): Quantitative solutions in hydrogeology and ground- water modeling.–Lewis Publishers, New York, 461 p. LAROCQUE, M., MANGIN, A., RAZACK, M. & BANTON, O. (1997): Contribution of correlation and spectral analyses to the regional stu dy of a large karst aquifer (Charente, France).– J. Hydrol., 205, 217–231. doi: 10.1016/S0022­1694(97)00155­8. LAROCQUE, M., BANTON, O. & RAZACK, M. (2000): Transient Sta­ te History Matching of a Karst Aquifer Ground Water Flow Model.– Ground Water, 38/6, 939–946. LUKAČ REBERSKI, J. (2008): Hidrogeološka i hidrogeokemijska os- nova za defi niranje slijeva Gacke i zaštita njenog izvorišta [Hydro geological and hydrogeochemical basis for the catchment area defi ni­ tion: River Gacka springs catchment area – in Croatian].– Unpubl. Master’s Thesis. Rudarsko­geološko­naftni fakultet, Sveuči li šte u Zagrebu (Faculty of mining, geology and petroleum en gi neer ing, University of Zagreb), 118 p. LUKAČ REBERSKI, J., KAPELJ, S. & TERZIĆ, J. (2009): An estima- tion of groundwater type and origin of the complex karst catchment using hydrological and hydrogeochemical parameters: A case study of the river Gacka springs.– Geol. Croat., 62/3, 157–178, doi: 10.4154/ gc.2009.15. LUKAČ REBERSKI, J. (2011): Određivanje podsljevova izvorišta rijeke Gacke na osnovi hidrogeoloških parametara [Definition of the river Gacka springs subcatchment areas on the basis of hydrogeological parametars – in Croatian].– Unpubl. Doctoral Thesis. Rudarsko­ geološko­naftni fakultet, Sveučilište u Zagrebu (Faculty of mining, geology and petroleum engineering, University of Zagreb), 158 p. MANDIĆ, M., BOJIĆ, D., ROLLER­LUTZ, Z., LUTZ, H. & BRONIĆ KRAJCAR, I. (2008): Note on the spring region of Gacka River (Croatia).– Isotop. Envir. Health Studies, 44/2, 201–208. doi:10.1080/ 10256010802066364. MANGIN, A. (1981): Utilisation des analyses corré1atoire et spectrale dans l’approche des systémes hydrologiques.– Computers Rendues Acad. Sci. Paris, 293, 401–404. Lukač Reberski et al.: Definition of the river Gacka springs subcatchment areas on the basis of hydrogeological parameters Geologia Croatica 53 MANGIN, A. (1984): Pour une meilleure connaissance des systemes hy- drologiques a partir des analyses correlatoire et spectrale.– J. Hy- drol. 67, 25–43. OFFICIAL GAZETTE NO. 91 (2008): Strategija upravljanja vodama [Water Management Strategy – in Croatian]. OFFICIAL GAZETTE NO. 47 (2008): Pravilnik o zdravstvenoj isprav­ nosti vode za piće [Drinking water regulations – in Croatian]. PADILLA, A., PULIDO­BOSCH, A. & MANGIN, A. (1994 Relative importance of Basefl and Quickfl ow from Hydrographs of Karst Spring.– Ground Water, 32/2, 267–277. PANAGOPOULOS, G. & LAMBRAKIS, N. (2006): The contribution of time series analysis to the study of the hydrodynamic character- istics of the karst systems: Application on two typical karst aquifers of Greece (Trifilia, Almyros Crete).– J. Hydrol., 329/3–4, 368–376. doi: 10.1016/j.jhydrol.2006.02.023. PARLOV, J. (2010): Identifikacija parametara za modeliranje toka pod­ zemne vode glavnih izvora u porječju Mirne [Identification of Pa- rameters Used for Groundwater Modeling of Major Springs in the Mirna River Valley – in Croatian].– Unpubl. Doctoral Thesis. Ru­ dar sko­geološko­naftni fakultet, Sveučilište u Zagrebu (Faculty of mining, geology and petroleum engineering, University of Zagreb), 111 p. PAVIČIĆ, A., KAPELJ, S. & LUKAČ, J. (2003): The influence of the Highway on the protected spring of Gacka river.– RMZ – Materials and geoenvironment, 50/1, 289–292. PLUMMER, L.N., PRESTEMON, E.C. & PARKHURST, D.L. (1994): An interactive code (NETPATH) for modelling net geochemical reac tions alog a flow path, Version 2.0, USGS.– Water­Resources Investigations Report 94–4169, Reston, Virginia. POSAVEC, K., BAČANI, A. & NAKIĆ, Z. (2006): A Visual Basic Spread­ sheet Macro for Recession Curve Analysis.– Ground Water, 44/5, 764–767. doi: 10.1111/j.1745­6584.2006.00226.x. POSAVEC, K., PARLOV, J. & NAKIĆ, Z. (2010): Fully Automated Objec tive-Based Method for Master Recession Curve Separation.– Ground Water, 48/4, 598–603. doi: 10.1111/j.1745­6584.2009. 00669.x. PRELOGOVIĆ, E. (1989): Neotektonski pokreti u području sjevernog Velebita i dijela Like [Neotectonic movements in northern Velebit and part of Lika – in Croatian]. Geol. vjesnik, 42, 133–147. RIGGS, H.C. (1968): Techniques of Water-Resources Investigations of the United States Geological Survey, Chapter A1, Some statistical tools in hydrology. Book 4, 39 p. STROJ, A. (2010): Podzemni tokovi u zaleđu krških priobalnih izvora na području Velebitskog kanala [Underground water flows in the hinterland of the velebit channel coastal karst springs – in Croatian]. Unpubl. Doctoral Thesis. Rudarsko­geološko­naftni fakultet, Sve­ učilište u Zagrebu (Faculty of mining, geology and petroleum en- gineering, University of Zagreb), 259 p. TALLAKSEN, L.M. (1995): A review of baseflow recession analysis.– J. Hydrol., 165/1–4, 349–370. doi: 10.1016/0022­1694(95)92779­D. TERZIĆ, J., STROJ, A. & FRANGEN, T. (2011): Hydrogeologic inves- tigation of karst system properties by common use of diverse meth- ods: a case study of Lička Jesenica springs in Dinaric karst of Croa­ tia.– Hydrological processes. doi:10.1002/hyp.91 TOMASKO, D., FISHER, A., WILLIAMS, G. P. & PENTECOST, E.D. (2001): A statistical study of the hydrological character of the Ed- wards aquifer.– Argonne National Laboratory, University of Chicago for the U.S. Department of Energy, 34 p. doi:10.1061/40856(200)149. VANDENSCHRICKA, G., VAN WESEMAELA, B., FROTA, E., PU- LIDO­BOSCH, A., MOLINAB, L., STIEVENARDC, M. & SOU­ CHEZD, R. (2002): Using stable isotope analysis (δD–δ18O) to cha- racterise the regional hydrology of the Sierra de Gador, south east Spain. – J. Hydrol., 265, 43–55. doi: 10.1016/S0022­1694(02)00097­5. VLAHOVIĆ, I., TIŠLJAR, J. & VELIĆ, I. (2009): Tercijarne karbonatne breče (paleogen­neogen­ Pg, Ng) [Tertiary Carbonate Breccia, PG- Ng – in Croatian]. – In: VELIĆ, I. & VLAHOVIĆ, I. (2009): Tumač Geološke karte Republike Hrvatske 1:300000 [Explanatory Notes of the Geological Map of the Republic of Croatia in 1:300.000 Scale]. Croatian geological survey, 141 p. ZANINOVIĆ, K., SRNEC, L., PERČEC­TADIĆ, M. (2004): Digitalna godišnja temperaturna karta Hrvatske [A Digital Annual Tempera- ture Map of Croatia – in Croatian].– Hrvatski meteorološki časopis (Croatian Meteorological Journal), Zagreb, 39, 51–58. ŽUGAJ, R. (2000): Hidrologija [Hydrology – in Croatian].– Rudarsko- geološko­naftni fakultet, Sveučilište u Zagrebu (Faculty of mining, geology and petroleum engineering, University of Zagreb), 407 p. Manuscript received October 09, 2012 Revised manuscript accepted November 27, 2012 Available online February 28, 2013 Geologia Croatica 66/1Geologia Croatica 54