COBISS: 1.01 TIME SERIES ANALYSIS APPLICATION FOR KARST AQUIFER CHARACTERISATION IN PINDUL CAVE KARST SYSTEM, INDONESIA APLIKACIJA ANALIZ ČASOVNIH VRST ZA KARAKTERIZACIJO KRAŠKEGA VODONOSNIKA Z JAMSKIM SISTEMOM PINDUL, INDONEZIJA Afid NURKHOLIS1, Tjahyo Nugroho ADJI1*, Eko HARYONO1, Ahmad CAHYADI1 & Slamet SUPRAYOGI2 Abstract UDC 551.44:556.33(594) Afid Nurkholis, Tjahyo Nugroho Adji, Eko Haryono, Ahmad Cahyadi & Slamet Suprayogi: Time series analysis application for karst aquifer characterisation in Pindul Cave karst system, Indonesia The karst hydrologic system in Pindul Cave can be categorised as 'binary' because it is recharged by the developed cave systems around it and by the Kedungbuntung ponor, which originates from an undeveloped karst system. The main objective of this study is to determine whether there are any differences in the aquifer characteristics that recharge to the Kedungbuntung ponor and the Pindul Cave outlet, using a time-series analysis approach. This analysis was conducted using univariate (autocorrelation) and bivariate (cross-correlation) statistical methods, which analyse data based on the domains of time (time-based analysis) and frequency (frequency-based analysis). Furthermore, the Master Recession Curve (MRC) calculation was conducted to confirm the time-series analysis result. Water level and rainfall data in both locations were recorded every 15 minutes over a six-month period (January-June 2017; N = 17,376). According to the time-series analysis, both locations were recharged by conduit, fissure and diffuse flow components; this observation was also confirmed by the MRC graphs which bore three components of recession. Comparison between the two locations showed that the Pindul outlet releases storage aquifers faster than Kedungbuntung. This condition may be caused by Izvleček UDK 551.44:556.33(594) Afid Nurkholis, Tjahyo Nugroho Adji, Eko Haryono, Ahmad Cahyadi & Slamet Suprayogi: Aplikacija analiz časovnih vrst za karakterizacijo kraškega vodonosnika z jamskim sistemom Pindul, Indonezija Kraški hidrološki sistem jame Pindul lahko označimo kot binarni, saj ga napajajo okoliški razviti jamski sistemi in ponor Kedungbuntung, ki se steka z nerazvitega kraškega sistema. Glavni cilj te študije je z uporabo analiz časovnih vrst ugotoviti, ali obstajajo kakršne koli razlike v značilnostih vodonosnika, ki napaja ponor Kedungbuntung, in tistim, ki napaja jamo Pindul. Ta analiza je bila izvedena z uporabo univariatnih (avtokorelacijskih) in bivariatnih (navzkrižno korelacijskih) statističnih metod, ki analizirajo podatke na podlagi časovnih in frekvenčnih vrst. Za potrditev rezultatov analiz časovnih vrst smo izvedli tudi izračun glavnih recesijskih krivulj (GRK). Vsakih 15 minut v šestmesečnem obdobju (januar-junij 2017; N = 17.376) so bili na obeh lokacijah zabeleženi podatki o vodostajih in padavinah. Glede na analizo časovnih vrst obe lokaciji napajajo komponente kanalskega, razpoklinskega in difuznega toka; to opažanje so potrdili tudi grafi GRK, ki kažejo na tri recesijske komponente. Primerjava med obema lokacijama je pokazala, da jama Pindul hitreje sprosti uskladiščeno vodo kot Kedungbuntung, kar je posledica kompleksnega kraškega sistema v kombinaciji z napajanjem iz ponora Kedungbuntung. Hkrati pa medsebojna korelacija praznjenja za 1 Karst Research Group, Faculty of Geography, Universitas Gadjah Mada, Indonesia, e-mail: afidnurkholis@gmail.com, adji@geo.ugm.ac.id, e.haryono@geo.ugm.ac.id, ahmadcahyadi@geo.ugm.ac.id 2 Department of Environmental Geography, Faculty of Geography, Universitas Gadjah Mada, Indonesia, e-mail: ssuprayogi@ugm.ac.id * Corresponding author Received/Prejeto: 20.03.2018 DOI: https://doi.org/10.3986/ac.v48i1.6745 ACTA CARSOLOGICA 48/1, 69-84, POSTOJNA 2019 AFID NURKHOLIS, TJAHYO NUGROHO ADJI, EKO HARYONO, AHMAD CAHYADI & SLAMET SUPRAYOGI the Pindul Cave outlet having a complex karst system supplemented by a developing system, in combination with recharge from the Kedungbuntung ponor. Meanwhile, the discharge-discharge cross-correlation between the two locations was less able to describe the complexity of the karst aquifers, as can be demonstrated by the tracing test method. In addition, the MRC indicated that the aquifers in the Pindul Cave karst system (which have developed) may still have a large storage capacity. To obtain more comprehensive information, especially with regard to the character of the developed karst system around the Pindul karst outlet, a follow-up study must be carried out over a longer period using several additional methods to support the results of this time series analysis study. Key words: Master Recession Curve, time series analysis, karst aquifer characteristics, Pindul Cave karst system. obe območji manj učinkovito opiše zapletenost kraških vodon-osnikov, kar se lahko dokaže s sledilnim poskusom. Poleg tega je MRC pokazal, da ima vodonosnik kraškega sistema jame Pindul še vedno veliko zmogljivost skladiščenja. Za pridobitev celovitejših informacij, zlasti glede značilnosti razvitega kraškega sistema okrog jame Pindul, bi bilo v nadaljevanju treba izvesti dolgotrajnejšo študijo, v kateri bi z uporabo več dodatnih metod ovrednotili rezultate te študije z analizo časovnih vrst. Ključne besede: glavna recesijska krivulja, analiza časovnih vrst, značilnosti kraškega vodonosnika, jamski sistem Pindul. INTRODUCTION Karst aquifers have three types of voids, which makes them complicated. Unlike other aquifers, the karst aquifer is challenging for assessment because of its porosity, discontinuity and the high variability of its hydraulic parameters (Ford & Williams 2007). Stevanovic (2015) also notes that recognising and understanding a karst aquifer development system involves an extended and difficult process because of its complexity and heterogeneity. White (1969) and Atkinson (1977) explain that the karst aquifer has diffuse porosity (the laminar flow that passes through homogeneous pores) and conduit porosity (the turbulent flow that passes through developed voids). Moreover, Atkinson (1985) highlights a kind of fissure porosity, which flows through medium voids. Ford and Williams (2007) and Rashed (2012) explain that there are different methods to characterise the karst aquifer, such as water balance estimation, borehole analysis, spring hydrograph analysis, water tracing, hydrochemograph analysis and aquifer modeling. In addition to those mentioned above, another method of identifying karst aquifer characteristics is to analyse the karst spring's response to the recharge (precipitation or concentrated infiltration) that flows via the ponors into the aquifer system. This method is possible because the spring, as an outlet of the karst system, is a feature that encompasses the majority of karst aquifer characteristics and processes (Bonacci 1993; White 2003; Adji & Misqi 2009; Adji 2012; Adji et al. 2017). The karst aquifer response is the spring's discharge, while the recharge is the precipitation or concentrated infiltration via the ponors. Spring hydrograph analysis can be performed by observing the aquifer response to determine the outflow components that contribute to a karst aquifer system because precipitation is the primary recharge. Furthermore, precipitation can be correlated with spring discharge. Time-series correlation between precipitation and spring discharge can provide an overview of aquifer response to precipitation as a recharge or input into the system. This approach reveals the true nature of karst aquifers in a particular geographical area. According to Larocque et al. (1998), time-series analysis is easily implemented and is often used to describe the characteristics of a karst aquifer. Time series analysis was carried out using both univariate and bivariate statistical calculations; these are differentiated, according to the unit of measurement, as either time-based or frequency-based. Auto-correlation (univariate) and cross-correlation (bivariate) analyse the data with respect to the time domain and can be converted to the frequency domain. For the latter, auto-correlation was translated into spectral density. Moreover, cross-correlation was changed to cross-amplitude as well as phase, coherency and gain functions. Explanations of the concepts and calculations of each time series analysis method can be found in Jenkin and Watts (1968), Man-gin (1984) and Box et al. (1994). Time-series analysis of karst aquifers was first practiced by Mangin (1981; 1984) and Mangin and Pulido-Bosh (1983) and was, in those times, based on a method of single variation, characterised by the individual structure of a time-series, with the correlation analysis characterised by the transformation of the input (rainfall) to the output (spring discharge). Since then, time series analysis studies utilising various types of hydrological data have been carried out by Padilla and Pulido-Bosch (1995); Eisenlohr et al. (1997); Petric (2002); Rahnemaei et al. (2005); Panagopoulos and Lambrakis (2006); Valdes et 84 ACTA CARSOLOGICA 48/1 - 2019 TIME SERIES ANALYSIS APPLICATION FOR KARST AQUIFER CHARACTERISATION IN PINDUL CAVE KARST SYSTEM, INDONESIA al. (2006); Herman et al. (2009); Liu et al. (2011); Zhang et al. (2013), and Kovačič and Ravbar (2016). Time series analysis can also be applied to a karst binary system; namely, one that receives autogenic recharge from another karst system. In other words, 'the karst aquifer is a part of the karst system' (Bakalowicz, 2005). Larocque et al. (1998) and Kovačič (2010) used temperature and electrical conductivity (EC) data (other than precipitation, water level, and discharge) in their time series analysis of karst binary systems. Additionally, Bailly-Comte et al. (2008) examined flood interactions between surface rivers and subsurface flow using time series analysis, while Mayaud et al. (2014) applied the same type of analysis and a MODFLOW model to assess overflow in binary karst systems. Thus, the purpose of this study was to apply time-series analysis to determine the difference in aquifer recharge between the Kedungbuntung ponor and the Pindul Cave outlet. Both locations are a unified system that can be categorised as karst binary. However, Kovacic (2010) and Kovacic and Ravbar (2016) claim that the classification of karst aquifers based on time series analysis alone does not adequately describe actual conditions and that other analyses are required to strengthen the results. As such, MRC analysis was also conducted in this study, in accordance with a review by Katsanou et al. (2015). GEOLOGICAL AND GEOGRAPHICAL SETTING The Pindul Cave karst system is located in the central re- described by Danes (1910); Lehmann (1936); Haryono gion of the island of Java. According to Van Bemmelen (2000); and Haryono and Day (2004). This subzone con- (1970), this karst system is located in Wonosari Basin sists of Miocene limestone of the Wonosari and Kepek (Fig. 1), stretching north of the Gunungsewu karst, as formations. The Wonosari formation comprises enor- Fig. 1: Regional geology of the Pindul Cave karst system (after Van Bemmelen 1970 and Kusumayudha 2005). ACTA CARSOLOGICA 48/1 - 2019 71 AFID NURKHOLIS, TJAHYO NUGROHO ADJI, EKO HARYONO, AHMAD CAHYADI & SLAMET SUPRAYOGI mous areas of coral-reef limestone in the south and bedded chalky limestone in the north (Balazs 1968; Van Be-mmelen 1970; Waltham et al. 1983; Surono et al. 1992). On the other hand, the limestone of the Kepek formation is bedded and comprises claystone, sandy marl and cal-carenite (Kusumayudha 2005). According to Rahaning-mas (2013), the local stratigraphy of the research area is made up of boundstone lithofacies and grainstone. Boundstone is a deposited coral-reef limestone from the Lower-Middle Miocene and is interfingered with the aforementioned grainstone. Grainstone itself is a carbonate limestone with coarse-to-smooth sand, deposited during the Middle Miocene. The Pindul Cave karst system has a passage length of approximately 300 m and a catchment area of 15.44 km2. Agniy (2016) classified the boundstone-dominated zone as a 'developed' karst aquifer (1.75 km2) and the grain-stone-dominated zone as an 'undeveloped' karst aquifer (13.69 km2) (Fig. 2). A developed karst aquifer is characterised by the presence of a number of connected cave systems - in this case, as Asri Cave, Greng Cave, Emas Cave, Candi Cave, Suruh Cave and Sioyot Cave, which extend from east to west (Haryono 2014). This underground river system has an outlet in Pindul Cave (Fig. 3) and the springs found in this area have a discharge range of 2-77 L/s, while some of the other underground river 458500 480000 4Ê1500 463000 484500 458500mE 460000 4E1500 4E3000 464500 - River ■ Rain Gauge H—Oyo strike-slip fault 1 Boundstone ^ Pindul Karst System Springs A A Pindul normal fault H Grainstone - Contour (meter) 0 Cave 0 Water level logger ■ ■ ■ Anticline fold axis ^J Packstone „ T-** ^ H *•»«' NiJSrV, rjrn Fig. 2: The Pindul Cave karst system (after Rahaningmas 2013 and Agniy 2016). 72 ACTA CARSOLOGICA 48/1 - 2019 TIME SERIES ANALYSIS APPLICATION FOR KARST AQUIFER CHARACTERISATION IN PINDUL CAVE KARST SYSTEM, INDONESIA Tab. 1: Discharge characteristics of caves and springs in the Pindul Cave karst system in the hydrologicalyear 2017. Pindul Cave outlet system 30,681 3,456 3,550 Kedungbuntung ponor 11,215 1,575 1,750 Emas - Greng - Asri cave system - - 650 Suruh - Candi cave system - - 210 Sioyot1 - Suruh - Sioyot2 cave system - - 780 Ngancar spring - - 77 Jebul spring - - 22 Beji spring - - 4 Mudal spring - - 10 Kali Banteng spring - - 2 Location Type Discharge (Q) L/s Maximum Minimum Average systems that fill the Pindul Cave outlet have a discharge range of 210-780 L/s (Tab. 1). Meanwhile, an undeveloped karst aquifer is characterised by the existence of an epikarst spring and surface river, which flows into the outlet of the Pindul Cave through the Kedungbuntung ponor (Fig. 3). According to Fig. 2, there is a relationship between the Kedungbun-tung ponor and the Pindul Cave outlet, with the ponor feeding the Pindul Cave system so that this karst system may be categorised as binary (Tab. 1). The rainfall data for a ten-year period (2006-2015), gathered by the Karangmojo rainfall station, shows that Fig. 4: a) Yearly rainfall data (in mm/year). b) Monthly rainfall data (2006-2015) at the Karangmojo rainfall station (in mm/month). ACTA CARSOLOGICA 48/1 - 2019 71 AFID NURKHOLIS, TJAHYO NUGROHO ADJI, EKO HARYONO, AHMAD CAHYADI & SLAMET SUPRAYOGI Tab. 2: The stage-discharge rating curves for the Pindul Cave outlet and the Kedungbuntungponor. . .. Number of discharge ,. , .. „2 Location measurement Stage-discharge rating curve R2 Pindul Cave (outlet of the karst system) 29 y = 6.01x43 0.96 Kedungbuntung (ponor) 29 y = 2,645x095 0.95 the research area has a monsoon climate, with the greatest monthly rainfall occurring in December (296 mm) during the October-April rainy season, while the dry season occurs between May-September (Fig. 4). The average annual rainfall is 1,725 mm, with the lowest level recorded in 2007 (926 mm) and the highest in 2013 (2,240 mm). DISCHARGE AND RAINFALL DATA In this study, water level and rainfall were measured every 15 minutes (N = 17,376). HOBO U20L-02 waterlevel data loggers were installed in the Pindul Cave outlet and in the Kedungbuntung ponor while an automatic rain gauge was installed in the centre of the karst area. The data recordings at both sites were conducted over a six-month period, from January to June 2017, and in order to subsequently convert this data into a discharge data series, stage-discharge rating curves were built, based on twenty-nine sets of discharge readings taken in both locations for peak, average, and low water-level conditions. The discharge measurements used the area-velocity method, and the measurement of flow velocity used a current meter. To maintain accuracy, flow velocity measurements were taken in straight-stream segments that were expected to have low turbulence. At each point, three flow velocity readings were taken in order to determine the average flow velocity. Then, the discharge was calculated using the mid-section method. Rating curves showed a good correlation between discharge and water levels at the different stations (R2 between 0.96 and 0.95; Tab. 2). Following this, discharge hydrographs were generated for the duration of the study period (Fig. 6). TIME-SERIES ANALYSIS Time-series analysis was applied in this research using univariate (auto-correlation) and bivariate (cross-corre- Fig. 6: Hydrographs of a) the Pindul Cave outlet and b) the Kedungbuntung ponor for a six-month period (January-June 2017). 74 ACTA CARSOLOGICA 48/1 - 2019 TIME SERIES ANALYSIS APPLICATION FOR KARST AQUIFER CHARACTERISATION IN PINDUL CAVE KARST SYSTEM, INDONESIA lation) statistical measurements. Both of these methods analyse data based on time but can be converted to the frequency domain. Auto-correlation can be changed into spectral density; meanwhile, cross-correlation can be transformed into cross-amplitude as well as phase, coherency and gain functions. A detailed explanation of the concept of time series analysis can be found in Jenkins and Watts (1968), Mangin (1984) and Box et al. (1994). In this study, time series analysis was calculated using XLSTAT software. Auto-correlation only uses one signal in its analysis. C(k=0.2) is a correlogram showing memory effect, which is a karst aquifer's ability to store an input signal (Mangin 1984). Zhang et al. (2013) compared forms of correlograms with flood recession hydrographs and identified those structures in the correlogram showing the duration release of both the diffuse (long-term response) and fissure flow components (intermediate-term response). The duration release of the fissure flow components is obtained from the initial sharp decline in the correlogram. Spectral density (S(f)) is a modification of auto-correlation in frequency form. S(f) shows the data distribution at various frequencies, where the peaks of the graph describe the flood cycle or pattern in the karst system. Treg (time regulation) quantified by spectral density indicates the storage capacity of a karst system (Larocque et al. 1998). Generally, a developed karst aquifer is believed to have a low memory effect and small storage capacity. According to Kovacic (2010) and Kovacic and Ravbar (2016), however, this view is distorted because a developed karst aquifer may also have a high memory effect and large storage capacity. Cross-correlation Cxy(k) is the relationship between the input signal and the output signal and is a cross-correlogram, identified by a flood hydrograph, which has a rising limb, peak discharge and recession period (Panagopoulus & Lambarakis 2006). The highest value of Cxy(k) denotes the average reaction of the time lag (Tlag) in a karst system. According to Zhang et al. (2013), the initial sharp decline in a cross-correlogram shows the duration release of the fissure flow component. Cross-amplitude (Sxy(f)) shows the filtering characteristics of the input signals by the karst system; when the S value is close to zero, this means that the T, re- xy lag sponse of the system in responding to the input signal (Zhang et al. 2013). An undeveloped karst aquifer dom-inantly filters input signals, so that the Sxy graph is close to zero over a longer period of time (lower frequency) in comparison to the developed aquifer. Phase (0 ) shows the difference in frequency between the input and output signals and occurs because the input signal is filtered by the karst system. This discrepancy also produces a graph similar to that of the flood hydrograph series. The mean of the time lag (Tlag) can be calculated from that series (Panagopoulus & Lambarakis 2006). The coherence function (COxy(f)) illustrates the linearity of the input-output in the karst system - a system is considered linear when the COxy value is close to or equal to 1 (Larocque et al. 1998). Linear relationships show that input directly impacts output in a developed karst system. Finally, the gain functions (Gxy(f)) show the releasing duration of the diffuse (base) flow when Gxy>1; the fissure (intermediate) flow duration when 0.4< Gxy<1; and the conduit (quick) flow duration when Gxy<0.4 (Padilla & Pulido-Bosch 1995). In a binary karst system, cross-correlation analysis and coherence functions need to be carried out to determine the hydrogeological interconnection between the two systems. Larocque et al. (1998) and Kovacic (2010) used discharge, water level, temperature, and electrical conductivity (EC) data to carry out this analysis. In the study outlined in this paper, discharge data were used to identify the hydrogeological connection between the Kedungbuntung ponor and the Pindul Cave outlet. MASTER RECESSION CURVE (MRC) A Master Recession Curve (MRC) was needed at each site to determine the release characteristics of the flow components in a karst aquifer over the course of a year. Posavec et al. (2006) implemented a visual basic program to calculate the MRC by way of an automatic matching strip method. In their study, the MRC was constructed using two recession equation types: (1) simple exponential for laminar (diffuse) flow; and (2) turbulent linear for turbulent (conduit) flow. According to Malik and Vojt-kova (2010), flood events may consist of several varieties of laminar and turbulent flow. Qt = Q0e_at (1) Qt=Q0(l-|3t) (2) where Qt is discharged at the time (t), Q0 is the initial discharge, e is an irrational number (approximately 2.71828), a is a recession constant for the diffuse flow, and p is a recession constant for the conduit flow. In this study, the MRC calculations in Pindul cave outlet and Kedungbuntung ponor were aided by RC 4.0 software (Gregor & Malik 2010). ACTA CARSOLOGICA 48/1 - 2019 71 AFID NURKHOLIS, TJAHYO NUGROHO ADJI, EKO HARYONO, AHMAD CAHYADI & SLAMET SUPRAYOGI RESULTS AND DISCUSSION DISCHARGE-DISCHARGE CROSS-CORRELATION ANALYSIS AND COHERENCE FUNCTION The discharge-discharge cross-correlation between the Pindul Cave outlet and the Kedungbuntung ponor had values of Tlag= 0 and C = 0.86 (Fig. 7); nevertheless, high correlation values and fast time lags do not indicate that the two locations are directly connected (Kovacic 2010). This value suggests that the Pindul Cave outlet and the Kedungbuntung ponor produced the same hydrologi-cal responses and time series structure. The hydrographs presented in Figure 6 show that both locations produced identical flood responses. According to Nurkholis et al. (2017), both locations experienced nine prominent flood events during the period from January to June 2017. In addition, both exhibited the same conduit flow T (1 lag v hour 45 minutes) as shown in the cross-correlogram (Fig. 9). In addition, the average discharge-discharge coherence function shows a far from linear value of 0.664 (Fig. 8). Larocque et al. (1998) explained that COxy values that are close to linear (0.75 and 0.81) indicate that the karst system is developed, while COxy values of 0.65 and 0.7 are indicative of one that is less developed. Meanwhile, Kovacic (2010) claimed that COxy with a value of 0.87 shows a karst system with a conduit tunnel that develops while a less developed karst has a COxy value = 0.45. In this study, the Kedungbuntung ponor recharges the Pindul Cave through a conduit tunnel that has been highly developed at the Tanding Cave (Fig. 9). Values of COxy that are far from linear are probably caused by the presence of sump in both the Kedungbuntung ponor -the Tanding Cave system and Tanding Cave itself - and the outlet of the Pindul cave system. This condition was confirmed by the results of the tracing test between the Kedungbuntung ponor and the outlet of the Pindul Cave (Agniy et al. 2017). A Tracer Breakthrough Curve (BTC) graph shows four peaks on different days, which indicates a complex system within two locations (Fig. 8). Fig. 7: a) Cross-correlation function and b) coherency function of discharge-discharge of the outlet of the Pindul Cave and the Kedungbun-tungponor. Fig. 8: a) Tracer Breakthrough Curves (BTC) of the Kedungbuntung ponor versus the Pindul Cave outlet (Agniy et al. 2017), b) the situation in Tanding Cave (July 28th, 2017; Photo: A. Nurkholis). 84 ACTA CARSOLOGICA 48/1 - 2019 TIME SERIES ANALYSIS APPLICATION FOR KARST AQUIFER CHARACTERISATION IN PINDUL CAVE KARST SYSTEM, INDONESIA However, in this study, the results of the time series analysis were unable to explain the complexity of the Ke-dungbuntung ponor and the Pindul Cave outlet, without the addition of the analysis obtained using the water tracing method. Meanwhile, Larocque et al. (1998) and Kovacic (2010) used electrical conductivity (EC) data in the cross-correlation method, which could then easily be compared to the tracer tests. Unfortunately, this study did not produce EC time series data. KARST AQUIFER RESPONSES BY TIME-SERIES ANALYSIS In this study, time-series analysis was used to describe the karst aquifer response to rainfall. This response was subsequently classified into the conduit (quick-flow), fissure (intermediate-flow) and diffuse (base-flow) compo- nents of discharge. Tab. 3 shows the Pindul Cave outlet and Kedungbuntung ponor responses in releasing their flow components. Conduit Flow Release The response of the conduit flow is shown by its average time lag (T, ) and maximum duration release. The T, o v lag' lag was obtained by a calculation involving the elements of cross-correlation, phase function and cross-amplitude, while the maximum duration release of the conduit flow was obtained by a calculation of the gain function. The result of the cross-correlation analysis show that the flow discharge was directly proportional to the rainfall. This relationship was demonstrated by Cxy, with the highest value (0.406 in the outlet and 0.366 in the Kedungbuntung ponor) of the Tlag = 1 hour 45 minutes (Fig. 9). Ac- Tab. 3: Flow component release in the Pindul Cave outlet and Kedungbuntung ponor. Karst system response Method Indicator Kedungbuntung Outlet Cross-Correlation time of correlogram peak 1 hr 45 min 1 hr 45 min Conduit Tlag Phase Functions e (f) hg 2irf 1 hr 15 min 1 hr 15 min Cross-Amplitude [S (f)] tends to be zero 1 hr 46 min 1 hr 46 min max duration Gain Functions gf 0,4 30 min 17 min Auto-Correlation Initially sharp decline 3 hr 20 min 3 hr 20 min duration Fissure Cross-Correlation Initially sharp decline 4 hr 4 hr duration range Gain Functions 0.4< g (f) <1 30 min - 6 hr 40 min 17 min - 1 hr 51 min Auto-Correlation C(k)= 0,02 >50 hr 6 hr 40 min Diffuse duration Spectral Density Function T = - • S(f = 0) reg 2 166 hr 40 min 41 hr 40 min min duration Gain Functions gxy(f)> 1 6 hr 40 min 1 hr 51 min Fig. 9: a) Cross-correlation functions and b) phase functions of precipitation (input) and discharge (output) at the Pindul Cave outlet and Kedungbuntung ponor. ACTA CARSOLOGICA 48/1 - 2019 71 AFID NURKHOLIS, TJAHYO NUGROHO ADJI, EKO HARYONO, AHMAD CAHYADI & SLAMET SUPRAYOGI cording to Panagopoulos and Lambrakis (2006), a small C value indicates that signal precipitation has been significantly reduced by the karst system. From the two values above, it appears that the Pindul Cave outlet had a greater Cxy value (0.406) than that found in the Kedung-buntung ponor (0.366), indicating that the karst aquifers around Pindul Cave are more developed. Furthermore, the results of the phase functions show the difference in frequency between input (rainfall) and output (discharge) caused by input filtering in the karst system. The graphic, which bears similarities to the flood hydrograph series, explains that filtering. Figure 9 shows the conduit flow frequency as 0.9x10-4-62.4x10-4 and the karst system outlet with a Tlag of 1 hour 15 minutes; the same as that found in the Kedungbuntung ponor. Cross-amplitude (Sxy) shows the input signal filtering characteristic by the karst system before release of the output signal (discharge). According to Zhang et al. (2013), this value is the maximum Tlag of the conduit flow and is obtained when the Sxy value is close to zero (Fig. 9). An undeveloped karst aquifer will filter the signal at a high frequency (short period) so that the Sxy value will be close to zero at a low frequency. The result shows that the Pindul Cave outlet and Kedungbuntung ponor have Sxy values close to zero and on the same frequency of 94x10-4 (Fig. 10). Therefore, both locations have a Tlag = 1 hour 46 minutes, based on the cross-amplitude graph. According to Padilla and Pulido-Bosch (1995), gain functions (Gxy) show the maximum duration release of the conduit flow when the value of Gxy is <0.4. Figure 11 shows that the Pindul Cave outlet had a value of G <0.4 xy at a frequency of >600x10-4; thus, a duration of <17 minutes. Meanwhile, the Kedungbuntung ponor had a maximum duration release in conduit flow at a frequency of >330x10-4; thus, a duration of <30 minutes (Fig. 11). The results showed that the aquifer recharging the Pindul Cave outlet releases a faster conduit flow than the Ke-dungbuntung ponor. In summary, the results of cross-correlation, phase functions, and cross amplitude analysis show that both locations had a similar Tlag characteristic - i.e., a high Tlag value (<1 hour 46 minutes). For comparison, a time series analysis by Kovacic (2010) using hourly data gave a T, of 17-37 hours. Furthermore, the MRC calculations lag demonstrated that both locations have complex conduit flow releases. This condition is indicated by the presence of two types of turbulent flow models (green and red), both at the Pindul Cave outlet and the Kedungbuntung ponor (Fig. 12). Meanwhile, the gain function analysis showed that Fig. 10: Cross-amplitude functions of precipitation (input) and discharge (outputs) at a) the Pindul Cave outlet and b) the Kedungbuntung ponor. Fig. 11: Gain functions of precipitation as input and discharge as output at a) the Pindul Cave outlet and b) the Kedungbuntung ponor. 78 ACTA CARSOLOGICA 48/1 - 2019 TIME SERIES ANALYSIS APPLICATION FOR KARST AQUIFER CHARACTERISATION IN PINDUL CAVE KARST SYSTEM, INDONESIA Fig. 12: The Master Recession Curve (MRC) of a) the Pindul Cave outlet and b) the Kedungbun-tungponor. both locations had slightly different duration releases in conduit flow (<17 minutes in the outlet and <30 minutes in the Kedungbuntung ponor). This condition was also confirmed by the MRC graph in the Pindul Cave outlet, which displayed a secondary turbulent flow (green) and was steeper than that of the Kedungbuntung ponor. This rapid release took place because when precipitation occurs, the Pindul Cave outlet is dominantly recharged by the cave system and the Kedungbuntung ponor. Meanwhile, the ponor is recharged only by the autogenic river. Fissure Flow Release The fissure flow release component was analysed according to the average and range of its duration. The average duration release of the fissure flow component is obtained from the initially sharp decline in the auto-correlation and cross-correlation graphs (Zhang et al. 2013). The duration range is calculated using the gain functions, which is when the graph shows values of 0.4< Gxy <1 (Padilla & Pulido-Bosch 1995). Auto-correlation and cross-correlation analyses showed that the Pindul Cave outlet and Kedungbuntung ponor had similar duration averages in their fissure flow components: approximately 3 hours 20 minutes and 4 hours (Figs. 9 & 13). The gain functions showed that the Pindul Cave outlet had a value of 0.4< Gxy <1 in a frequency of 90x10-4 - 600x10-4; thus, a duration range of17 minutes - 1 hour 51 minutes. Meanwhile, the Kedungbuntung ponor had a frequency range of 25x10-4-330x10-4; thus, a duration ranges from 30 minutes - 6 hours 40 minutes (Fig. 11). ACTA CARSOLOGICA 48/1 - 2019 71 AFID NURKHOLIS, TJAHYO NUGROHO ADJI, EKO HARYONO, AHMAD CAHYADI & SLAMET SUPRAYOGI Auto-correlation and cross-correlation of the data showed that both locations had the same average release times in their fissure flows. In addition, the duration range of the fissure flow at the Pindul Cave outlet had a faster release than the Kedungbuntung Ponor. This was also confirmed by the MRC graph in Figure 12, where the primary turbulent flow (red) corresponds to the release of the fissure flow. Here, the Pindul Cave outlet had a steeper primary turbulent flow so that the release of the component flow fissure was faster than that found in the Kedungbuntung ponor. Diffuse Flow Release Aquifer responses to diffuse flow are shown in terms of average and minimum durations. The average duration of the diffuse flow was obtained by implementation of the auto-correlation and spectral density functions, while the minimum duration was determined using gain function calculations. The auto-correlation calculation identifies the duration release of the diffuse flow from the average of all recession events. This duration average is obtained when the C(k) value reaches 0.2 on the graphic. This is also termed the 'memory effect' of the karst system (Mangin 1984). The Pindul Cave outlet C(k) value reached 0.2 within 6 hours 40 minutes, while the Kedungbuntung ponor did not attain this value within the maximum time period, as can be seen in the graph. This shows that the outlet had a more rapid diffuse flow release than the Ke-dungbuntung ponor (Fig. 13). The spectral density (S(f)) results also revealed a faster duration release of diffuse flow in the outlet than in the ponor. This faster duration release was confirmed by the time regulation value (Treg), which produced a frequency of 4x10-4 in the Pindul Cave outlet whereas the Kedungbuntung ponor emitted a frequency of 1x10-4. When converted to a duration domain, the outlet had a duration release in diffuse flow of 41 hours 40 minutes, which was faster than that found in the Kedungbuntung ponor (166 hours 40 minutes), as shown in Figure 13. For diffuse flow, Padilla and Pulido-Bosch (1995) stated that the minimum duration of diffuse flow occurs when the value of G >1. The Pindul Cave outlet reached xy this value at a frequency of <90x10-4 (duration >1 hour 51 minutes), while the Kedungbuntung ponor did so at a frequency of <25x10-4 (duration >6 hours 40 minutes), as shown in Fig. 11. From the auto-correlation, spectral density and gain function results, it was apparent that the Pindul Cave outlet released diffuse flow more rapidly than the Kedung-buntung ponor. Overall, the Pindul Cave karst system (which encompasses both locations) had a fast diffuse release or low memory effect. For comparison, the calculations of memory effect and regulation time carried out by Kovacic (2010) using hourly data produced values of 292 - 635 hours and 1,365 hours. This fast release occurred because the Pindul Cave karst system gains a recharge from the underground river system (a developed karst) in the vicinity and also from the autogenic river (the Kedungbuntung ponor). Meanwhile, the MRC (blue - Fig. 11) graph at both locations shows that diffuse flow storage tends to be stable, as indicated by a small constant value of a (0.004). According to Maillet (1905), the e-a value in equation (1) is a recession constant; therefore, diffuse flow recession constants in both locations may be calculated using the equation 2.71828-0004 to produce a value of 0.996. The recession constant calculation carried out by Adji et al. (2016) and Adji and Bahtiar (2016) showed that the value of 0.996 reflected stable diffuse flow storage throughout the year. This was evidenced by the flow conditions in the Pindul Cave outlet and Kedungbuntung ponor, which flow continuously over the course of a year, with an average discharge of 3,973 L/s and 2,195 L/s in the rainy season (January-April 2017); and 3,715 L/s and 1,805 L/s at the beginning of the dry season/end of the rainy season (May-June 2017). This result was in accordance with Kovacic's statement (2010); that is, a developed karst system can also have a large storage capacity. Fig. 13: a) Auto-correlation functions and b) Spectral density functions of discharge at the Pindul Cave outlet and Kedungbuntung ponor. 80 ACTA CARSOLOGICA 48/1 - 2019 TIME SERIES ANALYSIS APPLICATION FOR KARST AQUIFER CHARACTERISATION IN PINDUL CAVE KARST SYSTEM, INDONESIA CONCLUSIONS This study applied time series analysis to a binary karst system between the Pindul Cave outlet and Kedungbun-tung Ponor. The discharge-discharge coherence function calculation between these two locations produced a value of 0.664 (which was far from linear = 1). This condition showed that there was a contribution from the developed karst aquifer around Pindul Cave that recharged the outlet to the overall binary karst system. Furthermore, time series analysis may be used to analyse the characteristics of karst aquifers in releasing conduit, fissure, and diffuse flow components. The results of this study indicated that, in general, the Pin-dul Cave outlet had a faster release of flow components than the Kedungbuntung ponor. This is because aquifers around Pindul Cave have a more complex (developed) karst system structure. The Tracer Breakthrough Curve also demonstrated that the Pindul Cave outlet is recharged by two main systems: the developed cave system and the Kedungbuntung ponor. Moreover, the latter is an autogenic river that is recharged by karst aquifers and these tend to be undeveloped. Time series analyses are statistically based, so there is a possibility that study results may not reflect the actual physical functioning of a particular system. Accord- ingly, in this study, the MRC was used to confirm the time series analysis calculation. The existence of diffuse flow components was evidenced by the recession curve in the linear flow model, fissure flow in the primary turbulent flow model, and conduit flow in the secondary turbulent flow model. The results of the MRC analysis demonstrated that the Pindul Cave karst system tended to have stable storage of diffuse flow components. As a side note, the data used in this study had a short duration of just six months. In addition, the application of data from a specific hydrological year greatly influences time series analysis calculations. When data are used that span a longer time period, the calculations are more likely to reflect actual hydrological conditions. Furthermore, the results of this study need to be confirmed by follow-up studies - in particular, by installing hydrological measuring tools at Sioyot cave or Ngancar spring, so that the characterisation of the karst aquifers in the Pindul Cave system become more comprehensive. It is hopeful that a comparison in the recharge flows from the Kedungbuntung ponor and the underground river system around the Pindul Cave outlet (Asri, Greng, Emas, Candi, Suruh, and Sioyot caves) can be made using time series analysis and other methods. ACKNOWLEDGMENTS This research was funded by Hibah Penelitian Unggulan Perguruan Tinggi, The Ministry of Research, Technology and Higher Education of the Republic of Indonesia (scheme No. 015/SP2H/DPRM/II/2016). The completion of this manuscript is also being financed by Hibah Tim Pasca Sarjana No. 1780/UN1/DITLIT/DIT-LIT/ LT/2018. The authors would like to thank Romza Fau-zan Agniy and Hendy Fatchurohman for gathering the field monitoring data. Thanks also go to the blind reviewers who have provided input for improving this manuscript. REFERENCES Adji, T.N. & M. Misqi, 2009: The distribution of flood hydrograph recession constant for characterization of karst spring and underground river flow components releasing within Gunung Sewu Karst Region, Indonesia.- Indonesian Journal of Geography, 42, 1, 35-46. DOI: https://doi.org/10.22146/ijg.1575 Adji, T.N., 2012: Wet season hydro chemistry of Bribin River in Gunung Sewu Karst, Indonesia.- Envi- ron Earth Sci, 76, 1563-1572. DOI: https://doi. org/10.1007/s12665-012-1599-x Adji, T.N., Haryono, E., Fatchurrohman, H. & R. Okta-ma, 2016: Diffuse flow characteristics and their relation to hydrochemistry conditions in the Petoyan Spring, Gunungsewu Karst, Java, Indonesia.- Geo-sciences Journal, 20, 3, 381-390. DOI: https://doi. org/10.1007/s12303-015-0048-8 ACTA CARSOLOGICA 48/1 - 2019 71 AFID NURKHOLIS, TJAHYO NUGROHO ADJI, EKO HARYONO, AHMAD CAHYADI & SLAMET SUPRAYOGI Adji, T.N. & I.Y. Bahtiar, 2016: Rainfall-discharge relationship and karst flow components analysis for karst aquifer characterization in Petoyan Spring, Java, Indonesia.- Environmental Earth Sciences, 75, 735-745. DOI: https://doi.org/10.1007/s12665-016-5553-1 Adji, T.N., Haryono, E., Fatchurohman, H. & R. Oktama, 2017: Spatial and temporal hydrochemistry variations of karst water in Gunung Sewu, Java, Indonesia.- Environmental Earth Sciences, 76, 709-725. DOI: https://doi.org/10.1007/s12665-017-7057-z Agniy, R.F., 2016: Kajian hidrogeologi karst sistem Gua Pindul, Kecamatan Karangmojo, Kabupaten Gu-nungkidul.- Undergraduated Thesis, Universitas Gadjah Mada, Faculty of Geography, p.p. 86. (in Indonesian with English abstract) Agniy, R. F., Cahyadi, A. & A. Nurkholis, 2017: Analisis karakteristik akuifer karst dengan uji perunu-tan dan pemetaan gua.- Proceeding of 2nd annual conference of groundwater expert association (PIT-PAAI), 13rd-15th September 2017, Yogyakarta, p.p. 54. DOI: http://doi.org/10.17605/OSF.IO/DFXJH (in Indonesian) Atkinson, T.C., 1977: Diffuse flow and conduit flow in limestone terrain in The Mendip Hills, Somerset (Great Britain).- Journal of Hydrology, 35, 95-110. DOI: https://doi.org/10.1016/0022-1694(77)90079-8 Atkinson, T.C., 1985: Present and future directions in karst hydrogeology.- Ann. Soc. Geoi. Belgique, 108, 293-296. Bakalowicz, M., 2005: Karst groundwater: A challenge for new resources.- Hydrogeology Journal, 13 (1), 148-160. DOI: https://doi.org/10.1007/s10040-004-0402-9 Balazs, D., 1968: Karst regions in Indonesia.- In: Lang, S. (ed.) Karszt-Es Barlangkutatas. Globus Nyomda, 5, pp. 3-62, Budapest. Bailly-Comte, V., Jourde, H., Roesch, A., Pistre, S. & C. Batiot-Guilhe, 2008: Time series analyses for Karst/ River interactions assessment: Case of the Coula-zou river (southern France).- Journal of Hydrology, 349 (1-2), 98-114. DOI: https://doi.org/10.1016/'. jhydrol.2007.10.028 Bonacci, O., 1993: Karst springs hydrographs as indicators of karst aquifers.- Hydrologi-cal Sciences, 38 (1), 51-62. DOI: https://doi. org/10.1080/02626669309492639 Box, G.E.P., Jenkins, G.M. & G.C. Reinsel, 1994: Time series analysis: forecasting and control, 3rd ed.- Prentice Hall, pp. 598, New Jersey. Danes, J.V., 1910: Die karstphanomene in Goenoeng Sewoe auf Java.- Tijdschrift der Koninklijke Ned- erlandsche Aardrijkskundig Genootschap, 27, 247260. Eisenlohr, L., Bouzelboudjen, M., Kiraly, L. & Y. Ross-ier, 1997: Numerical versus statistical modeling of natural response of karst hydrogeological system.-Journal of Hydrology, 202, 244-262. DOI: https:// doi.org/10.1016/s0022-1694(97)00069-3 Ford, D.C. & P.W. Williams, 2007: Karst Hydrogeology and Geomorphology.- John Wiley & Sons Ltd, pp. 562, Sussex. Gregor, M. & P. Malik, 2010: RC 4.0 User's Manual Hydro Office Sofware for Water Science. [Online] Available from: https://hydrooffice.org/ Downloads?Items=Manual [Accessed 20th February 2019] Haryono, E., 2000: Some properties of epikarst drainage system in Gunungkidul Regency, Yogyakarta, Indonesia.- The Indonesian Journal of Geography, 32 (79-80), 75-86. Haryono, E. & M. Day, 2004: Landform differentiation within the Gunung Kidul Kegelkarst, Java, Indonesia.- Journal of Cave and Karst Studies, 66 (2), 62-69. Haryono, E., 2014: Speleogenesis Gua Pindul dan seki-tarnya.- Research Report, Universitas Gadjah Mada, Faculty of Geography.- [Online] Available from: http://karst.geo.ugm.ac.id/?page=risdet&iddet=61 [Accessed 20th February 2019] Herman, E.K., Toran, L. & W.B. White, 2009: Quantifying the place of karst aquifers in the groundwater to surface water continuum: A time series analysis study of storm behavior in Pennsylvania water resources.- Journal of Hydrology, 376, 307-317. DOI: https://doi.org/10.1016/jjhydrol.2009.07.043 Jenkins, G.M. & D.G. Watts, 1968: Spectral analysis and its applications.- Holden Day, pp. 525, San Francisco. Kovacic, G., 2010: Hydrogeological study of the Malenscica Karst Spring (SW Slovenia) by means of a time series analysis.- Acta Carsologica, 39 (2). DOI: https://doi.org/10.3986/ac.v39i2.93 Katsanou, K., Lambrakis, N., Tayfur, G. & A. Baba, 2015: Describing the Karst Evolution by the Exploitation of Hydrologic Time-Series Data.- Water Resources Management, 29 (9), 3131-3147. DOI: https://doi. org/10.1007/s11269-015-0987-x Kovacic, G. & N. Ravbar, 2016: Characterization of selected karst springs in Slovenia by means of a time series analysis.- In: Stevanovic, Z. (ed.) Karst without boundaries. CRC, Taylor & Francis Group, Leiden, pp. 131-154. Kusumayudha, S.B., 2005: Hidrogeologi karst dan ge-ometrifraktal di daerah.- Adicita Karya Nusa, Yo- 84 ACTA CARSOLOGICA 48/1 - 2019 TIME SERIES ANALYSIS APPLICATION FOR KARST AQUIFER CHARACTERISATION IN PINDUL CAVE KARST SYSTEM, INDONESIA gjakarta, Indonesia, p.p. 171. (in Indonesian with English abstract) Larocque, M., Mangin, A., Razack, M. & O. Banton, 1998: Contribution of correlation and spectral analyses to the regional study of a large karst aquifer (Charente, France).- Journal of Hydrology, 205, 217-231. DOI: https://doi.org/10.1016/s0022-1694(97)00155-8 Lehmann, H., 1936: Morphologische Studien auf Java.- Geographis-che Abhandlungen, pp. 114, Stutgart. Liu, L., Chen, X., Xu, G. & L. Shu, 2011: Use of hydrologic time series data for identification of hydrodynamic function and behavior in a karstic water system in China.- Hydrogeology Journal, 19, 1577-1585. DOI: https://doi.org/10.1007/s10040-011-0774-6 Maillet, E., 1905: Essais d'hydraulique souterrain et fluviale.- Librairie Scientifique, pp. 218, Paris. Malik, P. & S. Vojtkova, 2012: Use of recession-curve analysis for estimation of karstification degree and its application in assessing overflow/underflow conditions in closely spaced karstic springs.-Environmental Earth Sciences, 65, 2245-2257. DOI: https://doi.org/10.1007/s12665-012-1596-0 Mangin, A., 1981: Utilisation des analyses corrèlatoire et spectrale dans l'approche des systèmes hydrologiques.- Computers Rendues Acad. Sci. Paris, 293, 401-404. Mangin, A. & A. Pulido-Bosch, 1983: Aplicación de los análisis de correlación y spectral en el estudio de los acuiferos kársticos.- Tecniterrae, 51, 53-66. Mangin, A., 1984: Pour une meilleure conaissance des syste'mes hydrologiques a' partir des analyses corre'latoire et spectrale.- Journal of Hydrology, 67, 25-43. DOI: https://doi.org/10.1016/0022-1694(84)90230-0 Mayaud, C., Wagner, T., Benischke, R. & S. Birk, 2014: Single event time series analysis in a binary karst catchment evaluated using a groundwater model (Lurbach system, Austria).- Journal of Hydrology, 511, 628-639. DOI: https://doi.org/10.1016/jj.jhy-drol.2014.02.024 Nurkholis, A., Cahyadi, A. & R.F. Agniy, 2017: Spatial and Temporal Flood Characterization of Pindul Karst System, Gunungkidul Regency.- IOP Conf. Series: Earth and Environmental Science, 145. DOI: https://doi.org/10.1088/1755-1315/145/1Z012071 Padilla, A. & A. Pulido-Bosch, 1995: Study of hydro-graphs of karstic aquifers by means of correlation and cross-spectral analysis.- Journal of Hydrology, 168, 73-89. DOI: https://doi.org/10.1016/0022-1694(94)02648-u Panagopoulos, G. & N. Lambrakis, 2006: The contribution of time series analysis to the study of the hy- drodynamic characteristics of the karst systems: application on two typical karst aquifers of Greece (Trifilia, Almyros Crete).- Journal of Hydrology, 329 (3-4), 368-376. DOI: https://doi.org/10.1016/j. jhydrol.2006.02.023 Petrič, M, 2002: Characteristics of recharge-discharge relations in karst aquifer.- Carsologica, Inštitut za raziskovanje krasa ZRC SAZU, pp. 154,.Postojna, Ljubljana. Posavec, K., Bacani, A. & Z. Nakic, 2006: A visual basic spreadsheet macro for recession curve analysis.- Groundwater, 44, 764-767. DOI: https://doi. org/10.1111/j.1745-6584.2006.00226.x Rahaningmas, F.R., 2013: Geologi dan hubungan antara fasies karbonat dan jenis porositas terhadap pem-bentukan Gua Pindul, Desa Bejiharjo dan seki-tarnya, Kecamatan Karangmojo, Kabupaten Gunungkidul, DIY.- Undergraduate Thesis, Institute Science & Technology of AKPRIND Yogjakarta, Faculty of Mineral Technology, p.p. 69. (in Indonesian with English abstract) Rahnemaei, M., Zare, M., Nematollahi, A.R. & H. Sedghi, 2005: Application of spectral analysis of daily water level and spring discharge hydrographs data for comparing physical characteristics of karstic aquifers.- Journal of Hydrology, 311, 106-116. DOI: https://doi.org/10.1016/jjhydrol.2005.01.011 Rashed, K.A., 2012: Assessing degree of karstification: a new method of classifying karst aquifers. Sixteenth International Water Technology Conference Proceedings, 15 April 2012, Istanbul, 121-129. Stevanovic, Z., (ed.) 2015: Karst aquifers-characterization and engineering (in Professional Practice in Earth Sciences).- Springer International Publishing, pp. 692. Surono, Toha, B., Sudarno, I. & S. Wiryosujono, 1992: Geologi lembar Surakarta-Giritontro, Jawa.- Pusat Penelitian dan Pengembangan Geologi, Bandung, Indonesia , p.1407/6-1408/3. (in Indonesian) Valdes, D., Dupont, J.P., Massei, N., Laignel, B. & J. Rodet, 2006: Investigation of karst hydrodynamics and organization using autocorrelation and T-AC curves.- Journal of Hydrology, 329, 432-443. DOI: https://doi.org/10.1016/jjhydrol.2006.02.030 Waltham, A.C., Smart, P.L., Friederich, H., Eavis, A.J. & T.C. Atkinson, 1983: The caves of Gunung Sewu, Java.- Cave Science, 10 (2), 55-96. White, W.B., 1969: Conceptual models for carbonate aquifer.- Groundwater, 7 (3). DOI: https://doi. org/10.1111/j.1745-6584.2012.00923.x White, W.B., 2003: Conceptual models for karstic aquifers.- Speleogenesis and Evolution of Karst Aquifers, 1, 2-8. ACTA CARSOLOGICA 48/1 - 2019 71 AFID NURKHOLIS, TJAHYO NUGROHO ADJI, EKO HARYONO, AHMAD CAHYADI & SLAMET SUPRAYOGI Van Bemmelen, R.W., 1970: The geology of Indonesia. General Geology of Indonesia and Adjacent Archipelagoes.- The Haque, Government Printing Office Zhang, Z., Chen, Xi., Chen, X. & P. Shi, P, 2013: Quantifying time lag of epikarst spring hydrograph response to rainfall using correlation and spectral analysis.- Hydrogeology Journal, 21, 1619-1631. DOI: https://doi.org/10.1007/s10040-013-1041-9 84 ACTA CARSOLOGICA 48/1 - 2019