A comparison of regularized, sharp boundary and tear zone inversions along an MT profile in Sabalan geothermal field, Iran

This paper investigates magnetotelluric (MT) data recorded along a profile in the Sabalan geothermal region, NW of Iran. To find the range of relevant models consistent with the data, this study employed the so-called regularized, tear zone, and sharp boundary inversions. This study could effectively derive three alternative classes of models. Although the models show stable common resistive and conductive features there are some inconsistent details. Unaltered surface rocks and porous Basalt exhibit a high resistive overburden underlain by relatively more conductive Paleozoic sediments. A common signature of hydrothermal systems appears, and resistivities increase beneath a highly conductive clay cap in deeper parts. An intriguing feature resolved in the smoothest inversion model is a second deep conductor of 30 Ωm resistivities at a depth of 3 km, extending close to the surface. It can be related to the hot, solidified volcanic intrusions, resemblingthe heat source in a geothermal system. This study applied the two other inversion approaches for further hypothesis tests. Although the tear zone inversion re-establish the deep conductor (with 38 Ωm resistivities at 3 km depth), it is absent in the sharp boundary inversion result. This study concludes that the second deep conductor has a limited structure resolution.


Introduction
The MT sounding is widely used for geophysical exploration in geothermal fields [1,2], and the MT data inversion mainly provides its interpretation. "Smoothest model" approach, which minimizes data misfit and model roughness simultaneously, is the most common method for MT data inversion [3,4,5]. Although the smallest roughness constraint employed in this algorithm relieves the non-uniqueness problem of MT data inversion and stabilizes the solution, it could not clearly image the shape and sharp boundary of geological units in practical situations. Some 2-D inversion strategies have been suggested to solve this problem. Smith et al., (1999) and deGroot-Hedlin and Constable (2004) suggested sharp boundary inversion where a different model parameterization that allows for blocky structures is used. The models consist of a few layers with laterally varying thicknesses and resistivities ( [6,7] obtain a different class of inversion models where strong resistivity contrast between tear zones is allowed. In this method the model space is divided into different regions; tear zones where the model norm is minimized independently ( [8,9]).
MT measurements have been carried out in Sabalan geothermal field, NW of Iran ( Figure  1) to investigate the sub-surface electrical conductivity structure. Data analysis, 2D inversion results, and possible causes of enhanced conductivities have been discussed in EDC [10]. These measurements' main result is an anomalously conductive layer underlain by moderate resistivities at deeper parts, consistent with the clay alteration mineral in a high temperature geothermal field. However, 3D inverse solution results recognize a deep conductive zone related to the heat source [11]. In this paper, applied smooth, tear zone, and sharp boundary inversions of MT data along a profile in this region to define the range of relevant models consistent with the measured data, which are physically and geologically reasonable.

T data analysis
MT is a geophysical exploration method that measures the natural time-varying electromagnetic (EM) fields at the earth's surface to determine the electrical resistivity distribution beneath the subsurface. Impedance tensors (equation 1) and tipper vectors (equation2) are MT transfer functions calculated from the ratio of electric (E) and magnetic (H) fields components for a frequency range that determine the inductive response of a conducting sub-surface to an external time-varying magnetic field.
(where x, y and z denote perpendicular axis in horizontal and vertical directions in aCartesian coordinate system).
In the first step of MT data interpretation, dimensionality, strike, and distortionanalysiswereperformed to determinewhether the measured data are EM inductive responses of a large scale 2-D regionalgeoelectric structure superimposed by some local 3-D inhomogeneities or if the subsurface structure is 3-D in a regionalscale.

Dimensionality and strike analysis
β-skew angles computed from phase relationships contained in impedance data [12] show low values (less than 5˚) for most periods of all sites (figure 2) and reveal that the data may be regarded generally as 2D.

Distortion analysis
Furthermore, a multi-site-multi-frequency approach of Mc-Neice and Jones [13] provides a method to separate and remove galvanic distortion parameters (twist and shear angles) contaminated impedance responses of the regional geoelectric structure. The method employs a least square approach to fit the measured data with a seven-parameter model

Results and discussion
In regularized inversion of MT data a penalty function (E(m)) composed of the data misfit (Φd) and the model roughness (Φm) weighted by the LaGrange multiplier (τ) is minimized [3]: By the use of tear zone inversion approach, the smoothing function(Φm) is cancelled out in some predefined regions of the model space [8,9]. Alternatively, in the sharp boundary inversion, a different model discritization approach is used and 2D earth model is presented by laterally variable layers whose resistivities are described at a set of nodes [6,7].

Regularized Inversion
We applied a non-linear conjugate gradient algorithm by Rodi and Mackie [3] to obtain a smooth inversion model from the data set by removing distortion effects from measured data. The algorithm employs a second-difference operator as a stabilizing functional to produce an optimal minimum structure model. A trade-off parameter (τ) relatively weights the stabilizing function concerning the data misfit in the inversion penalty function. Numerical studies conducted to find the optimal τ value show that a value of τ=5 results in minimum RMS misfit as well as minimum model roughness. Figure (5a) shows the smoothest model obtained from the inversion of apparent resistivities and phases of both TE and TM modes whose error floors were 20% and 5%, respectively. The starting model was a homogeneous half-space with 100 Ωmresistivities. The resistivity model shows a resistive top layer, being more elevated between sites SAB021 and SAB025 which reflects unaltered surface rocks and overlies a conductive layer representative of Smectite-Zeolite alteration. From this conductor's base, a more resistive layer, indicating the reservoir, extends down to 1000 m a.

Tear Zone Inversion
To examine this deep conductor's robustness, this studyconducted a second set of inversions, incorporating a tear zone with an imposed resistive layer (1000 Ωm), which extends from -250 m depth down the model space bottom. This tear zone region was free to evolve during the inversion process. The result (figure 5b) shows that this resistive feature's imposition will produce a higher level of data misfit. Furthermore, the inversion model exhibits a conductive feature deeper than -250 m, suggesting that the inversion tries to re-establish the deep conductor, representing the heat source in the geothermal system.

Sharp Boundary Inversion
The depth to the geothermal reservoir layer as well as its thickness is of particular concern and led us to execute the third set of inversions. A five-layer-over-basement model with interface depths of 2000, 400, -2500, -5000, -9800 m a.s.l was adopted as the inversion staring model. Resistivities of the starting layers and the basement were 70, 7, 20, 600, 200, and 100 Ωm, respectively. This study set the parameter nodes below the measurement stations. The final SBI model obtained an RMS misfit of 1.88 is shown in figure (5c). By comparison, the smooth inversion result shows more lateral variations in resistivities both in the surface and deeper parts of the model. The top of the caprock is coincidentally well resolved in both models. However, the top of the reservoir layer's geometry and its thickness vary significantly between figures (5a) and (5c).

Conclusion
Three different classes of resistivity models were obtained by applying different approaches in model parameterization and regularization. Classical 2D smooth inversion results in a resistivity cross-section where a deep conductive structure with 30 Ωm resistivity was resolved at 3 km depth, as the heat source in a geothermal region. This study applied different regularization methods in a second approach by allowing sharp transitions of resistivity across tear zone borders. The result produces consistent images of the heat source. In the third class of modeling, this study derives a sharp boundary inversion model. Although the positions of boundaries can be recovered, more precisely, more numerical studies are essential to describe inconsistencies between the results of SBI and the two other classes of models.