http://journal.uir.ac.id/index.php/JGEET E-ISSN : 2541-5794 P-ISSN : 2503-216X Journal of Geoscience, Engineering, Environment, and Technology Vol 7 No 1 2022 Alpine, F. et al./ JGEET Vol 7 No 1/2022 7 RESEARCH ARTICLE Identification of Geothermal System In “Diana” Area, Indonesia Based On Magnetotelluric Data Modelling Fajar Alpine 1, Y. Yatini*1, Iqbal Takodama2 1 Geophysics Engineering, Universitas Pembangunan Nasional Veteran, Condong Catur 55283, Yogyakarta Indonesia. 2PSDG, Soekarno-Hatta Street, No.444, Pasirluyu, Regol, Bandung, West Java 40254. * Corresponding author : jeng_tini@upnyk.ac.id Tel.: +85-74-319-8262; Received: Jul 31, 2021; Accepted: Mar 30, 2022. DOI: 10.25299/jgeet.2022.7.1.7448 Abstract These days, the number of geothermal explorations is being increased to obtain a greater new potential of geothermal energy. One of the methods that is often used is magnetotelluric (MT). By MT, the components of a geothermal system can be delineated based on the resistivity values. This research’s main purpose is MT data modelling in 1 D and 2 D to delineate the geothermal system in the research area. There are 18 point of soundings, with a distance of about 1 – 3 km for each point. Bostick Transformation is used in 1 D modelling while Non-Linear Conjugate Gradient inversion is used as 2 D modelling with L – curve analysis as a method to obtain an optimal value of regularization parameter. Based on the analysis of 1 and 2 D models, the caprock zone was identified with a resistivity value of < 50 Ωm at a depth of 500 m with a thickness of about 250 m. The reservoir zone was identified with a resistivity value range of (50 – 100) Ωm located at a depth of 1000 with a thickness of about 500 m. Also, fault structures have been identified at the center of the research area. The regularization parameter used for the 2 D modelling is 5, which has obtained RMS values of 2.25% and 2.21% for each line. Keywords: Geothermal System, Inversion, Magnetotelluric, Resistivity 1. Introduction Indonesia had an increase in population from year to year, which has made Indonesia the 4th country in the world with the greatest population. Along with the increase of the population, the energy needs increased too. Energy supply should at least be 1.25% - 1.30% greater than the population growth. However, the energy availability in Indonesia is not proportional with Indonesia’s population which keeps increasing. In the case of 2012, there was an increase of 6.20% in population, while only 3.15% of increase in total primary energy supply (TPES) (Purnomo, 2014). Indonesia is at the meeting of three plates which are the Australian Plate, Eurasian Plate, and Pacific Plate, and perhaps also the Phillipine Sea Plate at the North of Sulawesi (Katili, 1975). The meeting of those plates has a role in the making of the most complex geology setting in the world. Indonesia has 129 active volcanoes, which spread along Indonesia’s islands. Geothermal energy has a connection with volcanic events, which is why lots of geothermal systems were created in Sumatra, Java, Southeast Nusa, Sulawesi, and Maluku. Indonesia is a country with 29 GW of geothermal reserve. Geothermal energy is planned to be utilized in a greater number to become one of the renewable energies to contribute as the national energy, reaching 4% by 2025 (Purnomo, 2014). One of the geophysical methods used for geothermal explorations is MT. This method used the propagation of electromagnetic wave with a great depth of penetration, which reached tens of metres to tens of kilometres (Vozoff, 1991). MT can represent the subsurface of a geothermal system based on the resistivity values. MT is effective to delineate the conductive layer compared to the resistive ones (Wameyo, 2005). This method is suitable for geothermal exploration because the caprock in a geothermal system can be delineated due to its conductivity which will have a contrast in values with the layer below it, reservoir. Before this research was conducted, an earlier research has been done in “Diana” Area by Kholid and Marpaung (2011). This research was conducted on 36 MT sounding points, which resulted in 2-D models. The results have identified the caprock zone with resistivity value of < 100 Ωm in the center of the area, which spreads to the North. The caprock zone is located in the elevation of – 1000 asl with a thickness of 500 – 1000 m. The reservoir zone is located below the caprock zone, which is at the depth of 1000 – 1500 m, which deepens to the East and West direction. The reservoir zone is bordered by a fault structure with East – West direction. In this research, Bostick Transformation will be used for 1 D modelling and Non-Linear Conjugate Gradient inversion for 2 D modelling. The results will be analysed based on resistivity values to identify the geothermal system components in Diana “Area”, which are caprock, reservoir, and fault, based on the contrast in the resistivity values shown. Then, the regularization parameter tau (τ) which gives the most optimal 2-D model will be decided. Last, both models will be compared and the most optimal model will be decided. 2. Geology and Stratigraphy http://journal.uir.ac.id/index.php/JGEET 8 Alpine, F. et al./ JGEET Vol 7 No 1/2022 The geothermal system in research is associated with Quaternary volcanic rocks. This area is characterized with volcanic rocks which has a composition of andesitic until trachitic. The morphology in the area is dominated with steep and corrugated hills, where cone – like shapes were found in few places. This cone – like shapes were probably the centre of former eruption of young volcanics which unfolded near the manifestation area. The corrugated hills showed the erosion that occurred from the older volcanic rocks. Volcanic activities in the research have occurred from the Tertiary age which was predicted to be undersea volcanoes which developed to be early Quaternary terrestrial volcanoes. Most of the products of Tertiary volcanics composed of andesitic until trachitic have eroded and caused the disappearation of the eruption source’s traces, and also caused intensive splitting which resulted in good enough permeability for fluid to flow. The next geology process is orogenesis which caused uplifting and resulted in the creation of land. Along the orogenesis process, volcanic activities still went on and created a cone – like volcanic at the Southwest of manifestation, with products of andesitic lava and lava breccia. Cone – like volcanic in the research area is predicted to be the heat source which has leftover heat from the magma chamber. The research area was dominated by volcanic rocks that consist of lava flow which spread widely, as well as volcanic domes. There are seven fault structures that are developing. A dilational fault jog was predicted to has occurred at the intersection of these faults, which created a media for hydrothermal fluidflow to the surface. These structures’ movement patterns are two strike – slip faults and five normal faults respectively. A depression is identified by a cliff which bordered said depression that has made a curved to half – radial shape. From the morphology and structures’ patterns, this depression is predicted to be the result of collapse of foldings that occurred before. The geothermal area in this research area is predicted to be bordered by this depression, where a manifestation in shape of hot spring has appeared. The stratigrapy of the research area consists of Andesit Porphyr Formation (Tp), Undivided Volcanic Formation (Tvt), Andesitic Lava Formation (Tlt), and Basaltic Andesitic Lava Formation (Tld), as shown in Figure 1. Based on these formations, the research area is dominated by andesite, andesitic lava, and basaltic andesite. There are also sheer joints that have created minor faults. There is also andesitic lava dome in the center of the research area, where the alteration process resulted in creation of chlorite and clays in the area. These formations occurred in the Tertiary age, which made the geothermal system in this area classified as non – volcanic geothermal system. Fig 1. Geology map of research area “Diana” (modified from Kholid dan Marpaung, 2011). Fig 2. Stratigraphy in research area, consists of Andesit Porphyr Formation (Tp), Undivided Volcanic Formation (Tvt), Andesitic Lava Formation (Tlt), and Basaltic Andesitic Lava Formation (Tld) (modified from Kholid dan Marpaung, 2011). Alpine, F. et al./ JGEET Vol 7 No 1/2022 9 3. Methods MT method is a geophysical method that utilizes natural EM fields. This EM field is generated by various physical processes which are quite complex with a very wide frequency spectrum (10-5 -10-4) Hz. At fairly low frequencies (<1Hz), the solar wind containing electrically charged particles interacts with the earth's permanent magnetic field, causing variations in the EM field. Variations in the audio frequency range (>1Hz) are mainly due to meteorological activity in the form of lightning. Lightning that occurs somewhere causes EM waves that are trapped between the ionosphere and the earth (vawe guide) and propagate around the earth. (Vozoff, 1991). The dependence of electric-magnetic phenomena on electrical properties, especially the conductivity of the medium (earth) can be utilized for exploration purposes using the MT method. This is done by simultaneously measuring the variation of the magnetic field (E) and magnetic field (H) as a function of time. Information about the conductivity of the medium contained in the MT data can be obtained by solving Maxwell's equations using simple models. Maxwell's equations contain Faraday's, Ampere, Gauss's laws which are shown by equations (1) to (4) 𝛻 × �⃗� = − 𝜕�⃗� 𝜕𝑡 (Faraday law) (1) 𝛻 × �⃗⃗� = 𝐽 + 𝜕�⃗⃗� 𝜕𝑡 (Ampere law) (2) 𝛻. �⃗⃗� = 𝑞 (Gauss I law) (3) 𝛻. �⃗� = 0 (Gauss II law) (4) Wherein �⃗� : electric field (V/m), �⃗� : magnetic induction (W/m), �⃗⃗� : magnetik field (A/m), 𝐽 : current density (A/m2) dan 𝑞 : density of charge (C/m3) Further relationship between the displacement current �⃗⃗� , electric field �⃗� and magnetic field �⃗⃗� can be written eq (5) to (7). �⃗⃗� = 𝜀�⃗� (5) �⃗� = 𝜇�⃗⃗� (6) 𝐽 = 𝜎�⃗� (7) where ε : Electrical permittivity (F/m), μ : Magnetic permeability (H/m),σ : Conductivity of the medium (S/m) Maxwell's equations can be rewritten as eq (9)-(11). 𝛻 × �⃗� = −𝜇  �⃗⃗�  𝑡 (8) 𝛻 × �⃗⃗� = 𝜎�⃗� + 𝜀  �⃗�  𝑡 (9) 𝛻 ∙ �⃗� = 0 (10) 𝛻 ∙ �⃗⃗� = 0 (11) Research in “Diana” Area was done based on 18 points of MT soundings. The results of these MT soundings are MT time series data. The first process to be done is to convert the time series into frequency domain. SSMT2000 software is used for fourier transform and robust processing. Robust calculation processes Discrete Fourier Transform (DFT) into crosspower. This process also applied the robust calculation to reduce noise in the data. Before robust processing, the parameters used have to be set up which is the Robust Parameter file (PRM) (Rogers, 2005). The result of robust processing is pseudo – resistivity and phase curves in frequency domain for every frequency. The next step is crosspower selection using MTEditor software. Data read by this software is MT Plot file created by SSMT2000, which shows the pseudo – resistivity and phase curves, along with its crosspower for every frequency. Other than that, this software also shows the characteristics from every sounding data, such as pseudo – resistivity, phase, impedance, strike direction, coherency, and others, along the frequencies obtained. (Rogers, 2005). Crosspowers affected by noise are deleted automatically or manually. With that, the data quality can be improved. The auto – edit option in MTEditor will delete crosspowers that are too far from the average value. That is why it is recommended to begin the selection by auto – edit, then continued by editing manually. The main obyective in this process is to smooth the pseudo – resistivity and phase curves by deleting crosspowers that are affected by noise (Mwakirani, 2012). After that, the next smoothing process will use WinGLink software. There are three types of smoothing, being Sutarno one of them, which calculates smoothing of the pseudo – resistivity curve based on the phase curve. In general, Sutarno is used to ensure the consistency of pseudo – resistivity and phase values. Second is D+ Smoothing, which associates the smoothing calculation between resistivity and phase, so the solution of 1 D Earth which suits both parameters the most can be obtained. This calculation proves to be the most optimal for 2 D and 3 D data. Last is numerical smoothing, which calculates the resistivity and phase curves independently (Geosystem SRL, 2008). In this research, D+ Smoothing was used because of its calculation that associated the resistivity and phase calculation. The obtained smooth curves will be modelled in 1 D and 2 D. Bostick Transformation was used for 1 D modelling. Bostick’s concept is trial and error, where model will be fit based on the observed data. The results from all the 1 D models will be used to make a cross section. A cross section is where the resistivity values from all the 1 D models will be interpolated. Based on the cross sections for every line, the resistivity values between soundings will be obtained. Non-Linear Conjugate Gradient inversion was used for 2 D modelling. In 2 D modelling, regularization parameter that gives the most optimal model will be determined by L – curve analysing. The regularization parameter that will be determined is tau (τ). Based on 1 D and 2 D models, resistivity distribution maps in determined depths will be made. 4. Results and Discussion The result of 1 D modelling in Figure 3 shows that there are four layers under surface. Not all 1 D models have the same number of layers due to the heteroginity of resistivity and phase values. In general, for all soundings, the resistivity values increased along with greater depth. This is due to the more massive and compact the rocks are in greater depths. Rocks near the surface will have lower resistivity values because of its non – compact characteristic. This low compactness can be due to weathering process. Points of soundings that show the pattern of increasing resistivity along with depth are MT 10, MT 11, MT 12, MT 13, dan MT 14 for 1st line, and all the soundings except MT 15A and MT 16 for 2nd line. There are few 1 D models that do not show the pattern as explained before. In sounding point MT 8 and MT 9 in the 1st line, there is a decrease in resistivity values in the third 10 Alpine, F. et al./ JGEET Vol 7 No 1/2022 layer. This is due to the location of the soundings that are close to the hot spring manifestation in “Diana” Area. The decrease in resistivity values is due to fluid flow or temperature increase which resulted by the fluid flow around it. Then in MT 9A, it is shown that the resistivity value is still very low even though the depth has increased. This is due to the location of sounding which is really close to the manifestation area. This sounding point is predicted to be right above the fault that plays a role in the geothermal system in “Diana” Area as a media for fluid flow. Then in the 2nd line, resistivity values in MT 15A has a decrease in the third layer. Location of MT 15A is parallel to MT 8 from the 1st line, which is why the decrease of resistivity values is also due to fluid flow and temperature increase. MT 16 has a resistivity value that is still very low even though it has reached a greater depth, like MT 9A. So, it is predicted that there is a fault below MT 16 which has a role as a media for fluid flow. In 1 D and 2 D modelling, the value range used will be the same. The resistivity value ranges are divided into three, which are low, medium, and high. Table 1 shows the classification of these resistivity values, and the interpretation for each range. All interpretations done will be based on Table 1. Table 1. Classification of resistivity values in “Diana” Area No. Classification Resistivity Range (Ωm) Interpretation 1 Low < 50 Caprock zone 2 Medium 50 – 200 Reservoir zone 3 High > 200 Massive rocks Fig 3. 1 D model of sounding point MT 10. a) Apparent resistivity curve. Dots represent the observed data, line represent calculated data; b) Phase curve. Dots represent the observed data, line represent calculated data; 1 D model, where the vertical axis represent depth and horizontal axis represent resistivity. Fig 4. 1 D cross section of 1st line in “Diana” Area. 1st line has 8 sounding points with a distance of about 1 – 2 km. The result showed a) caprock; b) reservoir; c) fault and d) hot spring. Generally, the resistivity increased with depth. MT 9A showed low resistivity until a depth of 4 km, which is caused by the fault located under the hot spring. a b d c Alpine, F. et al./ JGEET Vol 7 No 1/2022 11 Fig 5. 1 D cross section of 2nd line in “Diana” Area. 2nd line has 10 sounding points with a distance of about 1 – 3 km. The result showed a) caprock; b) reservoir and c) fault. Generally, the resistivity increased with depth. The cross sections of 1 D models showed that in general, the resistivity values increased along with depth. Below sounding point MT 9A, there is low resistivity value of < 10 Ωm from the surface until a depth of 4 km. This is due to the location of manifestation below this sounding point. The low value pattern is predicted to be a fault which is a media for fluid flow in the geothermal system of “Diana” Area. This fault is shown by dotted line in Figure 4. The low resistivity value is associated by fluid with high salinity in a geothermal system (Llera, dll, 1990). However, fluid salinity value in geothermal system of “Diana” Area is unable to be determined due to limitations of supporting data. In 2D modelling, there is a process to determine the right regularization parameter by L – Curve analysing using Microsoft Excel. The L – curve analysing is done by doing 2D inversions multiple times with a different value of 𝜏 each time. From every result of those inversions, the value of roughness/RMS to 𝜏 will be plotted to find its maximum curvature. The location of the maximum curvature in the L – curve shows balance between the requirements of an optimum model, which has enough smoothness as well as low enough RMS value (Hansen, 1992). Fig 6. L – curve: Roughness/RMS versus parameter 𝜏. Value of 𝜏 = 5, is the 𝜏 chosen. Figure 6, show that maximum curvature occurred at 𝜏 = 5. So, the 2 D modelling done in this research uses parameter 𝜏 = 5. Figure 7 shows the 2 D model for 1st line. The inversion was done using TE and TM data, and 𝜏 of 5, as explained before. The amount of iteration used is 30. The inversion resulted in RMS error of 2.25 %. Low resistivity value in sounding point MT 10 until MT 14 near the surface is interpreted as caprock zone. According to Ussher (2000), a zone that is located at the top of a geothermal system has a low resistivity due to temperature factor (70 – 200°C) and abundance of conductive clay. Clay is the result of hyrothermal alteration that occurred in said temperature. The low value pattern below MT 11 looks declining compared with its surrounding. This is predicted to be the result from a collapse by an intrusion below it. Based on the geology data, the geothermal system in “Diana” Area was made in a depression zone, characterized by half radial cliff as shown in the model, which is predicted to be the results from foldings that have occurred before. Below sounding point MT 11, it is predicted that there is a fault resulted by said collapse, and it can be identified by the resistivity value contrast. The caprock zone is located at the centre of the line, going to Northeast direction. The zone with a medium range of resistivity is right below the caprock zone, and is predicted to be the reservoir zone. This zone only reaches a depth of about 1 km. At this zone, the temperature has increased and became higher compared to the caprock zone above it. Components of geothermal system that have higher temperature is characterize by higher resistivity compared to the conductive zone above it. This higher resistivity value is due to rocks’ matrix that is less conductive compared to the fluid inside it. At this zone, mineralization is dominated by results of low conductivity alteration. High temperature alteration can increase rocks’ resistivity value because of the changing process from smectite clay to illitic or chloritic clay (Ussher, etc. 2000). Other than that, porosity also decreases along with greater depth, and can result in higher resistivity. Based on geology data and early researches, this high resistivity value is interpreted to be andesite rocks which are still very compact. Then, at a depth of below 2 km, the increase of resistivity until 1000 Ωm has appeared. This increase is predicted to be caused by metamorphic rocks which are characterized to be more resistive. 2 D modelling for 2nd line uses the same parameters as the 1st line. Figure 8 is the result of 2 D modelling for 2nd line. This model has a RMS error of 2.21 %. In this model, it is shown that the interpreted parts have the same location as the model for 1st line, which is at the Northeast of the a b c c 12 Alpine, F. et al./ JGEET Vol 7 No 1/2022 line. There are four faults interpreted which are below sounding point MT 17A, MT 19, MT 20, and MT 21, identified by resistivity value contrast. This fault is predicted to be the result of uplift by an intrusion below, and caused collapse near the surface. This collapse is shown by the radial shape of low resistivity near the surface, which keeps declining at the centre of the line. Fig 7. 2 D model of 1st line in “Diana” Area. 1st line has 8 sounding points with a distance of about 1 – 2 km. The result showed a) low resistivity: caprock; b) medium resistivity: reservoir; c) fault and d) hot spring. Generally, the resistivity increased with depth. MT 9A showed low resistivity until a depth of 4 km, which is caused by the fault located under the hot spring. The high resistivity at 4 km depth at the Northeast is interpreted as compact andesite rocks. Fig 8. 2 D model of 2nd line in “Diana” Area. 2nd line has 10 sounding points with a distance of about 1 – 3 km. The result showed a) low resistivity: caprock ; b) medium resistivity: reservoir and c) fault. Generally, the resistivity increased with depth. The hi gh resistivity at 4 km depth at the Northeast is interpreted as compact andesite rocks. Radial shaped caprock zone and reservoir zone suited the geology data which stated that there is a leftover cliff that makes a radial shape (Kholid dan Marpaung, 2011). This radial shape is predicted to occur due to an intrusion which caused uplift to the layer above it. Below the caprock zone, is a zone with higher resistivity, reaching 150 Ωm. This zone has a low temperature. According to Ussher (2000), high temperature located at the top of a geothermal system is caused by low fluid saturation, little hydrothermal alteration, and a decrease in resistivity caused by temperature. Temperature effect in caprock zone has been explained by a research by Kristmannsdóttir (1979). In this research, it has been explained that at a temperature of 50 – 100°C, alteration process in geothermal system produces smectite and zeolite as alterated minerals that dominate and caused a rock to be conductive. When depth increases, resistivity increases along with it. This is due to at a temperature of 220 – 240°C, smectite and zeolite is replaced by chlorite as dominating alterated mineral. When temperature increases, epidote mineral dominate and caused resistivity to increase even more (Kristmannsdóttir, 1979). That is why, it can be seen in Figure 7 and 8 that a reservoir zone has higher resistivity compared to the caprock zone above it. Even though reservoir zone contains fluid which can decrease resistivity value, there is no mineralization which a b c d a b c c c c Alpine, F. et al./ JGEET Vol 7 No 1/2022 13 produces conductive minerals that occurs in the reservoir zone. Conduction by minerals will give lower resistivity compared to rocks’ pores that contain fluid (Flóvenz, 2005). High resistivity value in a depth of 1 – 5 km is caused by massive rocks which has low porosity. Other than rocks’ compactness, high resistivity is also caused by temperature. At high temperature of a geothermals system, non conductive minerals dominate more and caused resistivity to increase. This caused conduction by rocks’ surface and fluid in pores dominate more compared to conduction by minerals. Transition from smectite to chlorite – epidote zone occurs at about 230°C (Flóvenz, 2005). Alteration process depends more on temperature, but porosity and permeability also affect in this process. That is why, it can be predicted that massive rocks that have low porosity at depth of 1 – 5 km and high temperature have low permeability, and caused decrease in intensity of alteration. These massive rocks are interpreted as andesite rocks based on the geology data. In Figure 9, it is shown that low resistivity zone is located at the East. This shows that caprcock zone reaches a depth of 500 until 1000 m. At the West, there is also a low resistivity zone and is predicted to be the caprock zone as well. However, the caprock zone at the West only reaches a depth of about 500 m. The reservoir zone is shown by the medium resistivity zone that appeared at the same location of the caprock zone. The reservoir zone is predicted to reach a depth of about 1500 – 2000 m. At greater depth, there are massive rocks that caused high resistivity to appear. At the South, it is seen that there is a pattern that elongated with a direction of Northwest – Southeast. This pattern is interpreted as a fault that plays a role in the geothermal system of “Diana” field. This fault is identified by the discontinuity of medium resistivity at the depth of 500 – 2000 m, where the boundary of this medium resistivity makes a straight pattern wih Northwest – Southeast direction. The low resistivity zone that is interpreted as the caprock zone is located at the East (Figure 10). This caprock zone is at the depth of 500 until 750 m. The medium resistivity zone which is at the same location of the caprock zone at 1000 m is the reservoir zone. This reservoir zone only reaches a depth of less than 1500 m. At a depth of > 1500 m, high resistivity spreads more along with greater depth. This is caused by massiveness of rocks which keeps increasing with greater depth. At 2000 m, it can be seen that a medium resistivity value appeared, which keep decreasing even more until 3000 m at the Southwest. This is caused by the manifestation in the Southwest of “Diana” Area. This has also been explained at the discussion for 2 D model, where the low resistivity value is due to fluid flow. The decrease in resistivity is also caused by the increase in depth. At 3000 m, there is a discontiunuity between high resistivity and low resistivity which makes a straight pattern with Northeast – Southeast direction. This pattern is predicted to be the fault that plays a role in the geothermal system. The difference between the maps based on 1 D and 2 D can be seen from the value distribution. For the 2 D maps, it can be seen that the distribution is more directed to the line’s direction while for the 1 D maps, there are resistivity zones that gathered in a point or few points. This is because of the 1 D models that have different pattern value for every sounding point. The 2 D modelled laterally and resulted in lateral distribution on the maps. Based on the models and resistivity maps, this research stated that the 2 D model is more optimal to delineate the research target. The reason for this statement is because of the 2 D model is able to show the geology condition of the area. Also, the 2 D model also shows the target better laterally and more supporting in interpretation because there are no values that gathered in one point and caused closure like the 1 D maps. Fig 9. Resistivity distribution maps in depth based on 1 D model. A) 500 meter ; B) 750 meter ; C) 1000 meter ; D) 1500 meter ; E) 2000 meter ; F) 2500 meter dan G) 3000 meter. Fig 10. Resistivity distribution maps in depth based on 2 D model. A) 500 meter ; B) 750 meter ; C) 1000 meter ; D) 1500 meter ; E) 2000 meter ; F) 2500 meter dan G) 3000 meter 5. Conclusions The result of this research showed that the resistivity value obtained in “Diana” Area from MT data processing is divided into three ranges, low, medium, and high. The caprock zone with < 50 Ωm is located at depth of (500 – 750) m. The reservoir zone with (50 – 100) Ωm reaches depth of less than 1500 m. A fault is identified at the centre of the research area based on the value contrast. The high resistivity that reaches 1000 Ωm at below 2 km is interpreted as andesite rocks which are still very compact. Second, the regularization parameter (𝜏) used in 2 D modelling is 5, and is the most optimal 𝜏 based on the L – curve analysis. This 𝜏 obtained a RMS error of 2.25 % for the 14 Alpine, F. et al./ JGEET Vol 7 No 1/2022 1st line and 2.21% for the 2nd line. Last, the 2 D model is shown to be more optimal to delineate the geothermal system in “Diana” Area because it was able to represent the geology of the area and showed the lateral distribution of the resistivity value. Also, it do not show closures caused by gathered values in few points. Acknowledgments We would like to thank Pusat Sumber Daya Mineral Batubara dan Panasbumi (PSDMBP) for their cooperation in providing the data for this research. Also, we would like to thank Ms. Indriati Retno Palupi for her discussion throughout the process of this research. References Flóvenz, Ó.G., Spangenberg, E., Kulenkampff, J., Árnason, K., Karlsdóttir, R. dan Huenges E., 2005. The role of electrical conduction in geothermal exploration. Proceedings World Geothermal Congress 2005, Antalya, Turkey. Geosystem SRL. 2008. A Guide to using WinGLink. Milan: Geosystem SRL. Grandis, Hendra. 2009. Pengantar Pemodelan Inversi Geofisika. HAGI, Jakarta. Hansen, P. C., 1992. Analysis of discrete ill – posed problems by means of the L – curve, SIAM Rev, 34, 561 – 580. Jones, A., G. 1983. On the Equivalence of the “Niblett” and “Bostick” Transformations in the Magnetotelluric Method. J Geophys, 53, 72 – 73. Katili, J. A. 1975. Volcanism and plate tectonics in the Indonesian Island Arcs. Tectonophysics. 26, 165 – 188. Kauffman, A. A., dan Keller, G., V. 1981. The Magnetotelluric Sounding Method. Elsevier Scientific Publishing Company. Kholid, M. dan Marpaung, H. 2011. Survei Magnetotellurik Daerah Panas Bumi Lili-Sepporaki, Kabupaten Polewali Mandar, Provinsi Sulawesi Barat. Prosiding Hasil Kegiatan Pusat Sumber Daya Geologi. Kristmannsdóttir, H., 1979. Alteration of basaltic rocks by hydrothermal activity at 100 - 300°C. Developments in Sedimentology, 27, 359 – 367. Llera, F., J., Sato, M., Nakatsuka, K. dan Yokoyama, H. 1990. Temperature dependence of the electrical resistivity of water saturated rocks. Geophysics, 56, 576 – 585. Lichoro, C. M. 2015. Comparison of 1-D, 2-D and 3-D Inversion Approaches of Interpreting Electromagnetic Data of Silali Geothermal Area. Australia: Proceedings World Geothermal Congress 2015. Mwakirani, Raymond. 2012. Magnetotelluric (MT) Data Processing. United Nations University Geothermal Training Programme. Pedersen, J., dan Hermance, J., F. 1986. Least squares inversion of one – dimensional magnetotelluric data: An assessment of procedures employed by Brown University. Surveys in Geophysics, 8 (2), 187 – 231. Peri, V. P., dll. 2014. Shallow Geophysical Evaluation of the Transition Zone between the Guarani and Yrenda – Toba – Tarijeno Aquifer Systems (Argentine Gran Chaco). Revista Mexicana de Ciencias Geologicas, 31 (1), 76 – 92. Postendorfer, G. 1975. Principles of Magnetotelluric and Prospecting. Geopublication Associates, 1 (5). Purnomo, 2014. Ketahanan Energi Nasional 2014. Jakarta: Dewan Energi Nasional. Rodi, W., dan Mackie, R., L. 2001. Nonlinear Conjugate Gradients Algorithm for 2-D Magnetotelluric Inversion. Geophysics, 66 (1), 174 – 187. Rodriguez, O. D., Enriquez, O. C., Fucugauchi, J. U., dan Arzate, J. A. 2001. Occam and Bostick 1-D Inversion of Magnetotelluric Soundings in the Chicxulub Impact Crater, Yucatán, Mexico. Geofisica Internacional, 40 (4), 271 – 283. Rogers, Stuart. 2005. Data Processing User Guide. Canada: Phoenix Geophysics Limited. Rogers, Stuart. 2010. V5 System 2000 MTU / MTU-A User Guide. Canada: Phoenix Geophysics Limited. Simpson, F. dan Bahr, K. 2005. Practical Magnetotellurics. Cambridge: Cambridge University Press. Sofyan, Y., Ehara, S. dan Daud, Y. 2009. Decades of Indonesian Geothermal Energy Growth. J. Geotherm, 31(3), 167 – 176. Unsworth, Martyn. 2009. Introduction to Electromagnetic exploration methods. Geophysics 223, University of Alberta. Unsworth, Martyn. 2016. Theory of electromagnetic (EM) field propagation in the Earth. Geophysics 424, University of Alberta. Ussher, G., Harvey, C., Johnstone, R. dan Anderson E. 2000. Understanding The Resistivities Observed in Geothermal Systems. Proceedings World Geothermal Congress 2000, 1915 – 1920. Vozoff, Keeva. 1991. The Magnetotelluric Method in Nabighian, M.N., Electromagnetic Methods in Applied Geophysics, Tulsa, Okla., Society of Exploration Geophysics, 37, pt.B, p.641-711. Vozoff, Keeva. 1980. Electromagnetic Methods in Applied Geophysics. Australia: D Reidel Publishing Company, Geophysical Surveys 4, 9 – 29. Ward, S., H. dan Hohmann, G., W. 1988. Electromagnetic Theory for Geophysical Applications. M. N. Nabighian, Ed., Electromagnetic Methods, Theory and Practice, Society of Exploration Geophysicist, 1, 131 – 311. Woldesemayet, B. 2013. Processing and 1D Inversion of Magnetotelluric Data from Dubti Area in Tendaho Geothermal Field, Ethiopia. United Nations University Geothermal Training Programme, 933 – 955. © 2022 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/