http://journal.uir.ac.id/index.php/JGEET E-ISSN : 2541-5794 P-ISSN : 2503-216X Journal of Geoscience, Engineering, Environment, and Technology Vol 08 No 02-2 2023 Special Edition Special Issue from “The 1st International Conference on Upstream Energy Technology and Digitalization (ICUPERTAIN) 2022” Dwikorianto, T.B. et al./ JGEET Vol 08 No 02-2 2023 Special Issue from The 1st International Conference on Upstream Energy Technology and Digitalization (ICUPERTAIN) 2022 67 RESEARCH ARTICLE Ambient Noise Data Processing to Obtain Group Velocity for Subsurface Structure Identification: Preliminary Research in Hululais Geothermal Field, Sumatra, Indonesia Tavip Dwikorianto1, Yunus Daud2,*, Agustya Adi Martha3, Aditya A Juanda4 1Doctoral Program in Geothermal, Department of Physics, Faculty of Mathematics and Natural Sciences, Universitas Indonesia, In donesia. 2Doctoral Program in Geothermal, Department of Physics, Faculty of Mathematics and Natural Sciences, Universitas Indonesia, Indonesia. 3National Research and Innovation Agency, Bogor 16911, Indonesia 4PT Pertamina Geothermal Energy, Jakarta 10110, Indonesia. * Corresponding author : ydaud@sci.ui.ac.id Tel.:+62-211-7867-222; fax: +62-211-7867-222 Received: May 20, 2023. Revised : May 31, 2023, Accepted: June 10, 2023, Published: 31 Juli 2023 DOI: 10.25299/jgeet.2023.8.02-2.13883 Abstract Hululais area lies in the pull-apart basins of the Ketaun Segment and Musi Segment fault as a part of the Sumatra Fault Zone (SFZ). The boundary normal faults of pull-apart basins play an important role as major discharge zones for geothermal fluid because the extensional stress is concentrated in the boundary normal faults. In order to identify the geothermal reservoir structure in Hululais Geo thermal Field (HGF), we introduce the local-scale study of the Rayleigh wave group velocity structure using ambient noise tomography (ANT). The ANT studies were collected using 18 seismometers inside 12 km2 area with a spacing of 125 – 500 meters, deployed across the fault structure for 1 month. More than two thousand Rayleigh Green’s Functions are extracted by cross-correlation at available station pairs. Using the estimated green function in this preliminary research, the group velocity as a function of the period can measure the dispers ion curve by using multiple filter technique (MFT) and fast marching surface tomography (FMST) scheme to obtain group velocity images. The tomography result as group velocity image shows the subsurface Rayleigh wave structure variation. The NW-SE main structure is reflected by the contrast velocity structure between the central part and the north eastern-south western sides. The central part shows the low periods which are associated with low wave velocity However the margin of the central part shows the high velocity in all periods. The ANT studies have been efficient in time and cost, however useful in subsurface structure interpretation in Hululais Geotherma l Field. Keywords: Pull-apart basin; Rayleigh wave, Green function, Ambient noise tomography, Cross-correlation, Group velocity 1. Introduction The ANT is the passive seismic method that uses the signal permanent source from surface waves which is considered to be related to the interaction of ocean swells with the seafloor near coastlines so the advantage is that it does not depend on the earth-quake occurrence and can be recorded at any time and any location (Yang Y et al., 2008). Classical reservoir imaging, particularly for oil and gas exploration, is based on active seismic surveys that rely on artificial impulsive sources (explosions and vibrator trucks). However, such operations have a high cost compared to the overall budget of a deep geothermal project and might affect the social acceptance of the project, especially in urban areas (Lehujeur M et al., 2015). The ANT studies have been efficient in time and cost because the data will be recorded in somedays even in somehours so it will reduce cost respectively. Ambient noise seismic is used for the regional study of tectonic and basin, volcano, and geothermal in many fields of the world. It is applied in Indonesia for the study of some basins in Central Java (Zulfakriza et al., 2012); Western East Java (Martha et al., 2015), East Java (Martha et al., 2016), East Java and Bali (Martha et al., 2017); Bandung – West Java (Wuryani et al., 2019); Sulawesi Island (Panshori et al.,2019). The study of the volcanic area in Merapi and surrounding volcanic centers (Koulakov et al., 2016); The Agung-Batur Volcano Complex, Bali (Zulfakriza et al., 2020). Ambient seismic noise is also applied in geothermal studies, such as in Jailolo (Vita et al., 2022), and Wayang Windu (Wahida et al., 2018). Theoretical and experimental research has shown that the cross-correlation of ambient noise records from two receivers provides an estimate of the empirical green’s function between the receivers/seismometer (Snieder, 2004). In the previous study, the average spacing distance between station pairs ranged from 1 to 25 km. In this study we tried to use a tighter spacing compared to previous studies, from 0.5 – 1 km to be 0.125 – 0,5 km, in order to obtain tomography images at shallower depths in the geothermal field. 2. Geological setting Hululais Geothermal Field (HGF) is situated in Rejang Lebong District, Bengkulu Province, Sumatera Island, http://journal.uir.ac.id/index.php/JGEET mailto:ydaud@sci.ui.ac.id 68 Dwikorianto, T.B. et al./ JGEET Vol 08 No 02-2 2023 Special Issue from The 1st International Conference on Upstream Energy Technology and Digitalization (ICUPERTAIN) 2022 Indonesia. This field is located about 150 km from Bengkulu City towards the Bukit Barisan mountain in the north. Hululais is one of thirteen geothermal systems situated in pull-apart basins on Sumatra Island which is the special tectonic setting where the major strike-slip Sumatran fault and volcanic zone coexist. Hululais area lies in the southwestern margin of the dextral strike-slip fault couples of the Ketaun Segment to the north and Musi Segment to the south as a part of the Sumatra Fault Zone (Fig. 1). These segment movements form a basin approximately 15 km in length and 5 km wide. Fig. 1. Tectonic setting of Hululais There are at least three main trending geological structures direction, specifically Northwest-Southeast (NW-SE), North East-South West (NE-SW) and North-South (N-S) patterns (Fig. 2). Lineament structure are indicated from satellite imagery. The NW-SE lineaments orientation are stepping Suban Agung Fault faults, Nibung Fault and Semelako Fault, interpreted as Musi Segment. The N-S lineament orientation are Cemeh Fault, Gregok Fault and Nusuk Fault. The NW-SE lineament orientation, i.e., Manghijau Fault. These faults are a minor fault of Musi Segment formed by N-S compression stress from subduction Australian-Pacific plate movement. Fig. 2. Geological Map of Hululais The reservoir zone in below of 1500 meter depth has high temperature reaches 300 oC and high permeability in the southern up flow zone. Meanwhile in the northern out flow zone has a temperature below 215 oC that encountered by wells whose the permeability is very low. By plotted in FS-HSH diagram, the gas sample that taken from the fumarole and the wells indicate that the reservoir is water/liquid dominated system. 3. Data and Method Seismic waveform data was recorded by seismograph stations which the stations were equipped with the seismographic sensors of Lennartz 3D Lite, the digitizer of Taurus, the time synchronizer of a GPS antenna, and the power supply of the Accu battery. The seismograph has 3 components, i.e. vertical component (Z) and the horizontal component (N-S and E-W) which the Z component is used to construct the green function. The tube concrete 40 cm in length was set in the hole of the observation station as a seismograph sitting for reducing the local noise from trees, humans or animals, etc. Noise recorded by these stations was used to image the group velocity under the research area. Fig. 3. The distribution of seismographic station locations in the research area The main part of the network was composed of 141 stations settled on the grid in about 12 km square coverage area with an average interstation distance of 500 meters and 125 meters in a certain area. The seismic waveform data was recorded by 15 portable and 3 stationary seismographs. The 3 static stations are set in the margin area in the northern (N004), the southern (N143), and the western (N042). The measurement in 15 stations was done in 1 – 4 days and then relocated to the next station. Fig. 4. Flow chart of the ANT methodology The initial stage of the ANT method is data acquisition by using a seismometer. Then the ambient noise data processing procedure divides into four principal phases as shown in Fig. 4 : (1) single station data preparation, (2) cross-correlation and temporal stacking, (3) measurement of dispersion curves, and (4) tomography (starting by resolution test, including error analysis and selection of the acceptable measurements). The first phase of data processing consists of preparing waveform data from each station individually in SAC data. After the preparation of the daily time series, the next step Dwikorianto, T.B. et al./ JGEET Vol 08 No 02-2 2023 Special Issue from The 1st International Conference on Upstream Energy Technology and Digitalization (ICUPERTAIN) 2022 69 in the data processing scheme (Phase 2) is cross-correlation and stacking. Cross-correlation is performed daily in the frequency domain. After the daily cross-correlations are returned to the time domain they are added to one another, or ‘stacked’ into weekly, monthly, yearly, etc. time series. After the daily cross-correlations have been computed and stacked, the resulting waveform is an estimated green function. Using the estimated the Green function, the group velocity as a period function can be measured using the multiple filter technique (MFT). The study obtains dispersion group velocities of the surface waves by applying multiple band-pass filters to the impulse response in the frequency domain using do_mft. Program in Seismology’s do_mft allows for the numbers of different input parameters. This study sets the units to counts, the type to Rayleigh and the filter parameter (α) 12.5 for all inter- station distances. Suitable filter parameters (α) are needed to produce clear information and help the balance of resolution trade-off between the time and frequency domains. The study performs multiple filter analysis on periods of 0.5–12.5 s because this includes the highest energy band or the highest spectral amplitudes. The dispersion curve gives various color information (blue- green-red) at each period and group velocity. Red color represents the highest energy and has the smallest degree of ambiguity than green and blue. Finally, using the fast marching surface tomography (FMST) scheme to obtain group velocity images. FMST implements forward and inverse problems using the fast marching method (FMM) (Rawlinson and Sambridge 2005) and subspace, respectively. The combination of the FMM for the calculation of the forward problem and subspace methods for inversion provide strong and stable condition at tomographic imaging. FMM (Rawlinson and Sambridge, 2004) is a grid-based eiconal equation used in the first step inversion to define the problem in the future. Meanwhile subspace method works by projecting squares approximation of Φ to an n-dimensional subspace of the model space. The disruption given by: (1) Where A is the projection matrix of M × n, G is the Fréchet derivative matrix, γ is the gradient vector (γ = ∂ Φ / ∂ m and γ = W - mγ). Search direction is given by aj in which the first search direction a1 is appropriate with the steepest ascent. Once the model is estimated disturbance, the study can update the current model and then retrace the propagation paths using Fast Marching Method (FMM) scheme. 3.1 Cross-correlation The cross-correlation of noise on the surface of the earth began to be developed by Shapiro and Campillo (2004) demonstrated empirical cross-correlation calculations of seismic noise wavefield recorded on two seismographs. The time derivative of the ambient noise cross-correlation Cij (1, 2, t) between two different seismograph stations, station 1 (located at r1 recording component i) and station 2 (located at r2 recording component j) to the time domain Green’s function (TDGF) Gij(r1; r2, t). (2) The TDGF Gij(r1; r2, t) relates a unit concentrated impulse displacement in the direction i at r1 to the displacement response in direction j at receiver r2. The TDGF which comes from random noises that propagate from stations 1 to 2 produces a positive correlation time delay t and the time-reversed TDGF produces a negative correlation time delay −t (cited from Martha, et al., 2017). The data processing consists of preparing waveform data from each station individually. The recording data from the seismograph was in minisheet and ASCII data and then converted into SAC (Seismic Analysis Code) data to get a smaller size of the data. Fig. 5. a. Cross-correlation processing schemes in this research, for A and B pair stations (Martha et al., 2017), b. Results of cross-correlation between station N026 and other stations For cross-correlation processing data, the data series reading was made in 23 segments per day which in each segment has a 1-hour time window. The sampling frequency was set to 100 samples per second in a window. The cross-correlation process was done in the first segment between 2 stations and continued to the 23rd segment then stack all of the segments in a day for getting the result of the cross-correlation of the pairing station in a day (Fig. 5a). The result of the cross-correlation of the pairing station in the whole day's survey was derived by stacking all cross- correlation of the first to the end of the day (Fig. 5b). We used a bandpass filter with a frequency limit between 1 to 1.8 Hz to extract all pairs of stations in the study area. 3.2 Dispersion measurement After the daily cross-correlations have been computed and stacked, the resulting waveform is an estimated green function. Using the estimated green function, the group velocity as a period function can be measured using the multiple filter technique (MFT). Fig. 6. a. Example of the cross-correlation waveform at pair stations. b. dispersion curve. The dispersion curve (Fig. 6b) is derived from a cross- correlation waveform (Fig. 6a) by conversion from period amplitude to period-velocity series where warm colours indicate the highest quality data with larger amplitudes. Dots indicate all possible dispersion curves. 70 Dwikorianto, T.B. et al./ JGEET Vol 08 No 02-2 2023 Special Issue from The 1st International Conference on Upstream Energy Technology and Digitalization (ICUPERTAIN) 2022 4. Results and Discussion 4.1 Checkerboard test Before performing the inversion with real data, we investigated the capability of the array to solve the velocity structure at different periods through the checkerboard test. The resolution test by using a checkerboard should be done in order to get a good resolution of the tomography image. Parameters resolution tests of grid size, smoothing, and dumping was done by trial and error. The grid size is made in 0.3 km to 0.6 km and the dumping and smoothing parameters are made in the range of 1-50. Fig. 7. Resolution test by using the grid size of 0.3 km Fig. 7 shows the resolution test with the grid size of 0.3 km, the dumping of 10, and the smoothing of 10 shows that the velocity model in all periods is not reconciled with the initial model. Fig. 8. Resolution test by using the grid size of 0.4 km Fig. 8 shows the resolution test with a grid size of 0.4 km, the dumping of 10, and the smoothing of 10. The result is good in all periods because they reconcile with the initial model. Fig. 9. Resolution test by using the grid size of 0.6 km Fig. 9 shows the resolution test with a grid size of 0.6 km, the dumping of 10, and the smoothing of 10. The velocity model in the eastern is better than in other areas due to the interstation distance in the eastern is denser (~125 meters) than in other areas (~250 meters). The result is good in periods of 0.1 to 1.9 seconds but not in periods of 2.5 seconds. From the three modeling above, the grid size of 0.4 km and using the dumping of 10, and smoothing of 10 are used for the tomography process in building the group velocity map in some periods. 4.2 Group velocity The two-dimensional distributions of group velocities were computed based on the iterative tomographic inversion of the group-velocity data for all available station pairs. After selecting the highest quality dispersion curve measurements, we compute group velocity maps from 0.1 to 4 s period. The two-dimensional distributions of group velocities were computed based on the iterative tomographic inversion of the group-velocity data for all available station pairs. Fig. 10. Group velocity map by using a grid size of 0.4 km; dumping and smoothing of 10 in certain periods Fig. 10 shows the tomography process in period 0.1 to 0.4 s which indicating the distribution of low velocity in almost cover in the whole area, however the velocity increase simultaneously with the period increasing (Fig. 10). It is interpreted as soft formation in the upper zone with a period of 0.1 – 0.4 s and harder formation below it with a period of 1.3 – 2.2 s. The soft formation is interpreted as caprock, and the harder formation is the reservoir. The period of 0.7 – 1.0 s is interpreted as a transition zone from caprock to reservoir. In some areas, the low-velocity anomaly in a high period (> 1.3 sec) is interpreted as a permeable zone. 5. Conclusion Seismic ambient noise tomography by cross-correlation of seismic ambient noise using seismographic data from several portable seismograph deployments describes the condition of the subsurface structure of Hululais area. From the dispersion curve measurement, the short inter-station distance will cover the shallow formation and the longer inter-station distance will cover the deeper formation. And from the group velocity identification, the low velocity anomalies in the shallow zone are correlated with caprock and the higher velocity anomalies in the deeper zone are correlated with the reservoir. However, the low anomalies in the deeper zone are interpreted as the permeable zone. Acknowledgment Grateful and appreciative to PGE’s Exploration and Exploitation Division for supporting the data, information, and discussion. Also, the authors thank PT Pertamina Geothermal Energy for permission to publish this study in Dwikorianto, T.B. et al./ JGEET Vol 08 No 02-2 2023 Special Issue from The 1st International Conference on Upstream Energy Technology and Digitalization (ICUPERTAIN) 2022 71 the International Conference Up Stream Energy Technology and Digitalization. References Koulakov, I., Maksotova, G., Jaxibulatov, K., Kasatkina, E., Shapiro, N.M., Luehr, B.G. et al., 2016. Structure of magma reservoirs beneath Merapi and surrounding volcanic centers of Central Java modeled from ambient noise tomography. Geochem. Geophys. Geosyst 17, 4195–4211. doi:10.1002/2016GC006442 Lehujeur, M., Vergne, J., Schmittbuhl, J., Maggi, A., 2015. Characterization of ambient seismic noise near a deep geothermal reservoir and implications for interferometric methods: a case study in northern Alsace, France. Geothermal Energy 3, 3. DOI 10.1186/s40517-014-0020-2 Martha, A.A., Widiyantoro, S., Cummins, P.R., Saygin, E., Masturyono., 2015. Upper crustal structure beneath East Java from ambient noise tomography: A preliminary result. AIP Conf Proc 1658, 030009. Doi: 10.1063/1.4915017 Martha, A.A., Widiyantoro, S., Cummins, P.R., Saygin, E., Masturyono., 2016. Investigation of upper crustal structure beneath eastern Java. AIP Conf Proc1730, 020011. Doi: 10.1063/1.4947379 Martha A.A, Cummins P, Saygin E, Widiyantoro S, Masturyono., 2017. Imaging of upper crustal structur beneath East Java–Bali, Indonesia with ambient noise tomography, Geosci. Lett., 2017) 4:14 DOI 10.1186/s40562-017-0080-9 Panshori, A., Martha, A.A., Maryanto, S., 2019. Imaging the Velocity Structure of Rayleigh Wave in Sulawesi Island Using Ambient Noise Tomography, IJASRE 5(1), 85095. Doi: 10.31695/IJASRE.2019.33072 Rawlinson, N., Sambridge, M., 2005. The Fast Marching Method: An Effective Tool for Tomographic Imaging and Tracking Multiple Phases in Complex Layered Media. Exploration Geophysics 36(4), 341–350. Doi: 10.1071/EG05341 Saphiro, N.M., Campillo, M., 2004. Emergence of broadband Rayleigh waves from correlations of the ambient seismic noise. Geophysical Research Letters 31, L07614. doi:10.1029/2004GL019491 Snieder, R., 2004. Extracting the Green’s function from the correlation of coda waves: A derivation based on stationary phase. Phy. Rev. E 69, 046610. 10.1103/PhysRevE.69.046610 Vita, A.N., Zulfakriza, Z., Martha, A.A., Rohadi, S., Heryandoko, N., Milkerreit, C., 2022. Preliminary Result of Rayleigh Wave Tomography beneathJailolo Volcanic Complex, North Moluccas, Indonesia using Ambient Noise. J. Phys.: Conf. Ser 2243, 012024. doi:10.1088/1742-6596/2243/1/012024 Wahida, A., Wijaya, H., Yudistira, T., Sule, M.R., 2018. Ambient Noise Tomography for Geothermal Exploration, a Case Study of WWs Geothermal Field. International Symposium on Earth Hazard and Disaster Mitigation (ISEDM). AIP Conf. Proc. 1987, 020101. Doi: 10.1063/1.5047386 Wuryani, S.D., Yudistira, T., Widiyantoro, S., 2019. Surface Wave Tomography Using Seismic Ambient Noise Data for Subsurface Imaging beneath Bandung Basin, West Java and Its Surrounding. IOP Conf. Ser.: Earth Environ. Sci. 318, 012032. doi: 10.1088/1755- 1315/318/1/012032 Yang, Y., Ritzwoller, M.H., 2008. Characteristics of ambient seismic noise as a source for surface wave tomography. Geochemistry, Geophysics, Geosystems 9(2), Q02008. doi:10.1029/2007GC001814 Zulfakriza, Z., Nugraha, A.D., Widiyantoro, S., Cummins, P., Sahara, D.P., Rosalia, S. et al., 2020. Tomographic Imaging of the Agung-Batur Volcano Complex, Bali, Indonesia. Front. Earth Sci 8, 43. Doi: 10.3389/feart.2020.00043 © 2023 Journal of Geoscience, Engineering, Environment and Technology. All rights reserved. This is an open access article distributed under the terms of the CC BY-SA License (http://creativecommons.org/licenses/by- sa/4.0/). http://creativecommons.org/licenses/by-sa/4.0/ http://creativecommons.org/licenses/by-sa/4.0/ http://creativecommons.org/licenses/by-sa/4.0/