Sabalan, with an elevation of 4,811 meters, is the second-highest mountain located in the Alborz mountain range in northwestern Iran. Although it is an inactive stratovolcano, the overall system of faults around the Sabalan peak can be considered as two common groups of linear and arcuate faults. Faults with the linear structure are mainly observed along with the NW-SE and N-S directions at the southern part of Sabalan. Arcuate faults around the Sabalan volcano are almost coincided with the crater of the volcano and have led to the penetration of volcanic masses. In this study, magnetotelluric data modeling was performed on one profile (“namely the MM profile”). The profile includes both types of fault systems in the region. By assessing and interpreting the resistivity section of the MM profile by 2D magnetotelluric data inversion, we have clearly detected the faults by features such as displacement and sudden changes in resistivity anomalies forming conductive zones with low resistivity anomaly along with the profile. The MM profile has almost continuous conductive anomaly (<30 Ω-m) between elevations of 1000 to 3000 m except in the eastern side, where a resistive block is found from the surface down to lower depths. Eight faults were detected across the MM profile. Our final represented geology section is consistent with actual observations and 2D magnetotelluric data inversion and prior studies.

Keywords: Magnetotelluric, Stratovolcano, fault detection, inversion, resistivity section. filter

Faults are geological phenomena that cause the displacement of anomalies and surface and subsurface layers and fractures. A resistivity section resulting from magnetotelluric data inversion can illustrate a correct image of the layers and subterranean abnormalities to depths of several kilometers in a region. In some cases, in these resistivity sections, because faults can be a source of water accumulation due to the gap between the two walls, they can have resistance anomalies along the fault line. On the other hand, changes and displacements in layers and abnormalities can help analyze and determine the position of faults in the region. The magnetotelluric method is a frequency electromagnetic (EM) tool that uses natural variations in the Earth's Electrical and magnetic field as an input source. Variations in the Earth's natural fields provide information and the ability to study the Earth's electric substructure at great depths. The extensive frequency range of EM signals eradicates the problematic presence of conductive overburden or sampling frequencies, allowing deep penetration. Natural MT signals arise from a variety of natural currents, including solar winds and thunderstorms. The frequency range of MT data can be from 40 kHz to less than 0.0001 Hz. The data is acquired in a passive mode applying a combination of electric sensors and induction coil magnetometers and can detect resistivity changes to great depths.¹

Researchers such as Cagniard² and Tikhonov³ developed the theory underlying the magnetotelluric method independently in the1950's. They have observed that the electric and magnetic fields associated with telluric currents that flow in the Earth due to variations in the Earth's natural electromagnetic field should relate to each other in a certain way depending on the electrical characteristics of the Earth. The ratio of the horizontal electric field to the orthogonal horizontal magnetic field gives the electromagnetic impedance. An essential advantage of the MT method is that it concurrently measures the electric and magnetic fields in two perpendicular directions. The electric sensors determine the electric field obtained from measuring the voltage difference between electrode pairs E_x and E_y. Induction coils measure the magnetic field components in three orthogonal directions. The ratio of the recorded electric and magnetic fields gives an estimate of the Earth's apparent resistivity at any given depth.

Faults can move layers to great depths in the earth's crust. The MT method can be a beneficial approach for fault detection. The total faults of the Sabalan region, in general, have expanded in two forms of linear and arcuate.⁴ The cracks and fractures caused by these faults have caused the region to be very prone to geothermal reserves in the country. So far, faults in the Sabalan area have not been investigated using magnetotelluric data. However, in Iran and other countries, the mechanism of many significant faults has been investigated using magnetotelluric data.^5–7 In Sabalan, the first magnitude data was acquired in 1998 containing 212 soundings,⁸ and the study of the Sabalan geothermal area has several years of history. In the following years, the acquisition of magnetotelluric data in this region with a frequency range of 0.003 Hz to 320 Hz was carried out in 13 stations in order to investigate more Sabalan geothermal field. In this paper, using these valuable MT data, processing and performing accurate interpretation, including observations in the site, we seek to alginate the faults in the southern part of the Sabalan region.

Geology of the area

Geological problems of the Sabalan area, including tectonic conditions, structural dynamics, and geothermal, are mainly related to Sabalan volcano activities. In terms of tectonics, the studied area is part of the Sabalan slopes, which is part of Azerbaijan's highland plateau. The tectonic activity of the area is affected by tectonic plates, where the Caspian, Eurasian and Arabic plates encounter. The lithology variation of outcrops in the study area is related to rock unit formation in different volcanic phases. In Sabalan geological maps, especially in the studied area in the Moeil Valley, no significant rock outcrops have been shown, and alluvial deposits are predominantly observed. Sabalan's early activities have begun from Eocene, but the main activity that started the Sabalan Mountain began in Pliocene and continued until after the last frost. This volcano's material has been derived from deep magma but is influenced by complex processes, which the most important of them are digesting and disintegrating two magmas and crystallization. The evolution of Sabalan magma has taken place over a long time.⁹

Based on existing rock-time units in the area, the foreign experts of the New Energy Organization of Iran have divided the Sabalan region into four main units,¹⁰ introducing four formations with local names. The situation and arrangement of these formations are presented in Table 1.¹¹ The description of the four formations is as follows:

1. Dizue: The main part of the Dizue Formation consists of discontinuous Alluvial deposits with variable thickness. The most developed parts of this formation are in the Moeil and Dizou Valleys. The thickness of this alluvial set varies from 10 to 150 meters at different locations.
2. Kasra: The Kasra Formation is mainly associated with outcrops of andesite-trachea volcanic lava, mostly located in the upper parts of Sabalan, formed after forming the crater. The dispersion of the Kasra formation rocks is in the southeast that forms the upper parts of the Sabalan Mountains.
3. Tous: Most of the Tous formation consists of lava flows of dacite-trachea, andesite-trachea, and Rhyolite, which are dispersed in the form of dwarf domes and larches. The outcrops of the Tous formation are located in the northern parts.
4. Valhazir: Valhazir formation is formed from trachea-andesitic, tuff, and cavernous slab lava flows. The rocks of this formation with a thickness of at least 1000 m include trachea-andesitic lava, tuff, and pitcher cuttings. One of the essential features of the Valhazir formation is the development of thermal alteration in trachea-andesitic rocks.

The fault system and linearity are studied by.⁴ In general, they have categorized and grouped the faults in the region into two groups of linear faults structures and arcuate fault structures: For linear fault structures, most of Sabalan linear structures have NW-SE and NE-SW orientation. A significant and long fault in the region with NW-SE orientation located in the west and northwest of Sabalan is the leading cause for forming the Moeil Valley (with the trend of NW-SE). The arcuate structure of this arcuate fault around the Sabalan volcano is more or less in line with the crater of the volcano and is visible at the site of a sudden change in the topographic slope. The formation of this fault is associated with volcanic activity and has led to the penetration of volcanic masses. The set of these faults are shown in Figure 1 with a simplified geological map of the region. In Figure 2, the curve for fault stretches is shown for the study area. The main stretches are toward NW and NE. These fault sets are derived from geological maps and digital maps of Sabalan's geothermal field.¹²

Figure 1 Location of MT stations and simplified geological map of the Sabalan region.

Figure 2 Curve for fault stretches in the study area.

The magnetotelluric method is a passive surface geophysical procedure which uses natural Earth electromagnetic fields to investigate the structure of substrate electrical resistance.² In this method, electromagnetic waves acquire information from depths of the earth according to their frequency and specific resistance of the layers. Using the impedance tensor, which relates the horizontal components of the magnetic field to the horizontal components of the electric field, one can extrapolate the characteristic specific resistance and impedance phase to interpret the magnetotelluric data. The tensor analysis provides an impedance analysis of the degree of heterogeneity in the environment and distribution of the MT field with heterogeneity, which allows the determination of geological structures and the calculation of their size and extension:¹³

$\left(\begin{array}{c}{E}_x\\ E_y\end{array}\right)=Z \left(\begin{array}{c}{H}_x\\ H_y\end{array}\right)Z=\left(\begin{array}{cc}{Z}_{XX}& Z_{XY}\\ Z_{YX}& Z_{YY}\end{array}\right)$ (1)

Where Z represents the apparent impedance at a specified measuring station. can be calculated using eq.2.

$Z_{ij}=\frac{E_i}{H_j}{e}^{i\phi }\to \text{ φ}={\text{Ψ}}_{{\text{E}}_{\text{i}}}-{\text{Ψ}}_{{\text{H}}_{\text{i}}}={\tan }^{-1}\frac{Im\left(Z_{ij}\right)}{Re\left(Z_{ij}\right)}$ (2)

Where and ,  are the electric and magnetic field’s phases, respectively, and the impedance phase is obtained by eq.3, within a homogeneous field.

$Z=\frac{\omega {\mu }_0}{k}=\frac{\omega {\mu }_0}{{\left(-i\omega \mu \sigma \right)}^{\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$2$}\right.}}=\sqrt{\frac{i\omega {\mu }_0}{\sigma }}\stackrel{ }{\Rightarrow }\phi =\frac{\pi }{4}$ (3)

This tensor provides information for conductivity, strike and the dimension of conductor structures. The time series measured in different frequency bands are converted to the frequency domain. The transverse power spectrum is calculated to estimate the impedance tensor as a function of frequency. The determinant of the impedance tensor, known as effective impedance, is defined by eq.4.¹⁴

$Z_{det}=\sqrt{Z_{XX}{Z}_{YY}-Z_{XY}{Z}_{YX}}$ (4)

By using the determinant of impedance tensor, the apparent velocity is calculated (see eq.5) where  and  are magnetic permeability of free space and angular frequency respectively.

${\rho }_{\alpha }=\frac{1}{\omega {\mu }_0}{\left|Z\right|}^2=\frac{1}{\omega {\mu }_0}{\left|\frac{E_x}{H_y}\right|}^2$ (5)

In a two-dimensional ground, electrical conductivity is constant along a horizontal direction (if it is specified; the same direction as the Strike in the geological plane). The Maxwell equations are partitioned into two independent polarization modes of TE (electric field parallel to the geological Strike) and TM (flow current perpendicular to geological Strike). On Earth, a hypothetical coordinate system is considered, in which the x-axis is parallel with the strike of the structure (the y-axis is along with the profile), and the profile is perpendicular to the z-axis and faces downward. In this situation, the components of the magnetic field and the electric field in TE mode are Hz, Hy, and Ex, respectively. In this mode, the electrical currents are parallel to the subsurface structure. The magnetic and electrical field components described in TM mode are Ey, Hz, and Hx, respectively. In a one-dimensional and two-dimensional fields where the field data components are acquired along, and perpendicular to the strike, the components of the primary diameter of the impedance tensor have zero values. Generally, the TE mode is useful for detecting deep, and nonconductor structures and TM modes are useful for detecting shallow and conductive structures.¹⁵ Diagnosis of TE and TM modes was performed using polarization diagrams in each sounding for this study.

The components of the  impedance (Z) tensor are obtained from the mixed ratios of orthogonal elements of horizontal electric and magnetic fields in the frequency domain. In our case, all the measurement stations are located over a line. The data only allows the application of a two-dimensional interpretation process which requires identification of TE and TM modes agreeing to electric and magnetic fields parallel to the geological strike, respectively. Because the geological strikes are not known in priorly, the elements of electromagnetic fields are measured in geomagnetic or optional directions, and the impedance tensor is rotated to main axes. The direction of strike often changes with depth in the field and, therefore, the rotational angle changes at each frequency. For two-dimensional structures, there are various conditions and consequently, many potential schemes to determine the rotational angle.¹ Here, we have decreased the diagonal elements of the impedance tensor. There are two potential strike directions for a specific frequency, and the interpreter identifies the TE and TM modes employing geological and geophysical information. In the MT method, the properties of the rocks resistivity tensor are determined by measuring the natural time of vertical components in magnetic and electrical fields. In these measurements, two components of the electric field (Ex and Ey) and three components of the magnetic field (Hx, Hy, and Hz) are measured. Two sets of non-polarizing electrodes measure the electric fields, and two sets of perpendicular inductive coils measure the magnetic field. The arrangement of electrodes and electromagnetic coils in the MT data acquisition is presented in Figure 4. The magnetotelluric data was acquired in the Sabalan region with MTU-5A instrument with a frequency range of 0.003 to 380 Hz. The location of MT stations and the designed profile MM is shown in Figure 4. The MT data were acquired using three sets of MTU-5A equipment, two roving sites within Mt. Sabalan and one remote reference site. The raw time series data were processed using the Phoenix Geophysics, Ltd. SSMT2000 software¹⁶ and the resulting EDI files were edited, analyzed and modeled using the Geosystem’s WinGlink software.¹⁷ Overall, the MT data obtained were of good quality down to 0.01 Hz.

Figure 3 Arrangement of electrodes and coils in MT data acquisition.

Figure 4 The topography map including the location of stations, Faults and the MM Profile in the region. The black line shows the MM profile.

Data processing consists of two general steps. In the first step, MT data is transmitted to the frequency domain. The second step involves calculating the components of the impedance tensor, using frequency spectra. The response function can be used to calculate models to limit the noise. The data was acquired in several stages to limit the model parameters of response functions. These operations include a study on the initial obtained values, sample diversity and signal-to-noise ratio.¹³ In an area consisting of homogeneous or horizontal layers without noise, the apparent resistivity is easily calculated by measuring the amplitude of the Ex and Ey at one time and frequency and multiplying it with the square ratio of the obtained values. This operation can be performed for all frequency portions in order to determine the curve values.¹⁴ For final editing and resistivity modelling in this research, the WinGlink¹⁷ software was used. The MM profile with an approximate length of 12 km is located east-west at the southern part of the Sabalan volcano. The modelling was carried out on 13 stations along with it. Figure 5 shows the data of four stations, 220, 238, 244, and 248 in the form of resistivity and phase-frequency for the frequency range. During the acquisition of MT data, at some stations, due to bad weather conditions (i.e. snow, wind) or situated along the steep and rocky flanks of the Sabalan mountain, some data were edited or deleted from the impedance tensor elements. After these considerations, the data was used for inversion.

Figure 5 The elements of the impedance tensor Zxy and Zyx for stations 220, 238, 244 and 248. In each subplot, the blue squares and the red circles stand for yx and xy respectively.

Based on¹⁸ for dimensional analysis of MT data from impedance tensor, the skew values were calculated for all stations at required frequencies. The result is shown in Figure 6. A maximum error of 5% was applied to the impedance values of the 13 MT stations. This quantity shows the asymmetry in the environment. The results show that for most frequencies, this value is below 0.2, indicating that those structures are mostly 1-dimensional and 2-dimensional. The inversion accuracy is checked in Figure 7 by comparing the measured pseudo sections (field data) of phase and resistivity with calculated pseudo sections (model response). High diversity indicates a flaw in inversion results. The pseudo sections derived from inversion of MT data along the MM profile are shown separately for each of the TE and TM modes. These sections represent a good match between the field data and the model response, indicating an accurate in inversion. Dimensionality analysis of the measured data, based on the phase-sensitive skew,¹⁹ reveals that skew values are less than 0.3 for most stations and periods and the whole data set may be regarded as 2D (Figure 6). The phase tensor scheme²⁰ was employed to provide estimates of regional geo-electric strike direction for all profiles evinces that a significant strike direction of 0˚ is quite reasonable for the study area. Accordingly, 2D inversions were applied to the un-rotated impedance data Figure 7. In Figure 8, the resistivity sections are derived from 2D inversion of MT for each mode of TE and TM. As mentioned, each of these modalities can be investigated for examining different structures. In this paper, using the combination of both TE and TM modes, the final resistivity section was modelled, examined and interpreted.

Figure 6 Skew Parameters for 13 MT Stations.

Figure 7 Pseudo resistivity and phase sections derived from inversion of MT data on the MM profile, including TM and TE modes. a) Observed pseudo resistivity for TM mode, b) Observed pseudo resistivity for TE mode, c) Calculated pseudo resistivity for TM mode, d) Calculated pseudo resistivity mode for TE, e) Observed phase for TM mode, f) Calculated phase for TM mode, g) Observed phase for TE mode, h) Calculated phase for TM mode.

Figure 8 Resistivity section derived from inversion of MT data on the MM profile for a) TM method and b) TE method.

The position of the primary faults FLM1, FLM2, FLM3, FLM4, FLM5, FLM6, FLM7, and FLM8 that interrupt the MM profile is shown on the topographic map of the area (Figure 4). Figure 9 shows the results of the 2D inversion of MT data by combining TE and TM modes along with the MM profile to a depth of 4.5 km. The main characteristic of this section is the continuity of an anomaly (conductor) with a resistivity of lesser than 30 ohms at elevation 1000 to 3000 meters above sea level, except in the eastern part of the block where a high-resistivity block is present from deep depth to lower depths. This conductive anomaly reflects the geothermal flow that extends from the southwest of Sabalan to the Moyel Valley in the west of Sabalan.²¹ There is a perfect match between the fault structures in the region with outcrops of the faults on land and the result of 2D MT data inversion. In the final electrical resistivity section (Figure 9) FLM5 and FLM1 faults located near stations 238 and 23, respectively, have been detected due to the presence of resistivity anomaly between 30 to 100 ohms. Surface water penetrations from the surface to deep depths along the strike of these faults has made them visible in the resistivity section. The same feature is visible for the FLM6 fault located near station 232 with low resistivity contours.

Figure 9 The final electrical resistivity (TE+TM) along with MM profile.

Moreover, the FLM6 fault is detectable on the resistivity section with displacement and transportation of the low resistivity anomaly. Other faults (between station 6 and 244) detected by displacement and in some areas by the transportation of the low resistivity anomaly from depths 3000 above sea level to the surface. Figure 10 shows the geology of the region for the MM profile. This section is achieved by integrating the geological data of the region with the resistivity section resulting from an inversion of the MT data. According to the geological information of the area, the geothermal current observed in Figure 8 with a low resistivity of 20 ohms from ground level up to a height of about 1000 m above sea level can be due to the trachea-andesite alteration flow that extends from the southern part of Mount Sabalan to the western and northern parts of Sabalan. On the other hand, in the resistivity section, in deeper depths, an intrusive mass with a high resistivity exists. This high resistive mass can also be an intrusive mass that is associated with volcanic activity in this area, such as a dome of Kasra²¹ which is shown in the geologic section as an intrusive multiphase complex. The geological formations of the area are also well-marked in the geology section.

Figure 10 The final geology section derived by the inversion of MT data and geological evidence for the MM profile.

Since Faults can move the geological layers to a great depth inside the earth's crust, the MT method can be beneficial for their investigation. In this research, the inversion of magnetotelluric data has well illustrated the faults of the southern volcano of Sabalan. Modelling along a magnetotelluric profile including 13 MT stations, the resulting resistivity from the inversion of these data was derived up to a rough depth of 4.5 km. The system of faults in the Sabalan region is generally arcuate and linear. These faults cut in MT profile in 8 areas. Investigation and interpretation of MT resistivity profile clearly show that some of the faults in this region, including FLM5, FLM1, and FLM6 are detected with linear resistivity features along with the fault line. This feature is detectable due to the interface of surface water along with fault fractures. Other faults in the resistivity section are detected by sudden resistivity changes which coincide with the displacements of layers. The features and structures detected by the MT profile have consistency with the geology and field observations.

The data was provided by Institute of Geophysics University of Tehran (IGUT). We are very thankful for all of the great supports provided by the authorities of IGUT and the GeoMine Company of Iran. We also wish to thank the two anonymous reviewers of this manuscript.

We have no conflicts of interest to disclose for the review process of our manuscript.

1. DePaolis RA, Janota CP, Frank T. Frequency Importance for Words, Sentences and Continuous Discourse. Journal of Speech, Language and Hearing Research (JSLHR). 39(4): 714–723.
2. Beck DL, Flexer C. Listening is Where Hearing Meets Brain in Children and Adults. Hearing Science. 2011;18(2): 30–35.
3. Pereira-Jorge MR, Andrade KC, Palhano-Fontes FX, et al. Anatomical and Functional MRI Changes after One Year of Auditory Rehabilitation with Hearing Aids. Neural Plasticity. 2018:13.
4. Glick HA, Sharma A. Cortical Neuroplasticity and Cognitive Function in Early-Stage, Mild-   Moderate Hearing Loss: Evidence of Neurocognitive Benefit from Hearing Aid Use. Front Neurosci. 2020;14:93.
5. Beck DL, Clark JL. Audition Matters More as Cognition Declines and Cognition Matters More as Audition Declines. Audiology Today. American Academy of Audiology. 2009.
6. Beck DL, Harvey M. Issues in Cognition, Audiology and Amplification. Hearing Review. 2021.
7. Amieva H, Ouvrard C, Giulioli C, et al. Self‐reported hearing loss, hearing aids, and cognitive decline in elderly adults: A 25‐year study. J Am Geriatr Soc. 2015;63(10):2099–104.
8. Amieva H, Ouvrard C. Does Treating Hearing Loss in Older Adults Improve Cognitive Outcomes? A Review. J Clin Med. 2020;9(3):805.
9. Seeto M, Tomlin D, Dillon H. The Relations Between Auditory Processing Scores and           Cognitive, Listening and Reading Ablities. Ear Hear. 2021.
10. Golub JS, Brickman AM, Ciarleglio AJ, Schupf N, Luchsinger JA. Association of subclinical hearing loss with cognitive performance. JAMA Otolaryngol Head Neck Surg. 2020;146(1):57–67.
11. Gaeta L, Azzarello J, Baldwin J, Ciro CA, Hudson MA, Johnson CE, John AB. Effect of reduced audibility on      Mini-Mental State Examination scores. J Am Acad Audiol. 2019;30(10):845–855.
12. Shinn-Cunningham BG, And Best V. Selective Attention in Normal and Impaired Hearing. Trends Amplif. 2008;12(4):283–299.
13. Alikcovic E, Ng EHN, Fiedler L, et al. Effects of hearing aid noise reduction on early and late cortical representations of competing talkers in noise. Currently in revision for Frontiers in Neuroscience. 2021.
14. O’Sullivan J, Herrero J, Smith E, et al. Hierarchical Encoding of Attended Auditory Objects in Multi–Talker Speech Perception. Neuron. 2019;104(6):1195–1209.
15. Hausfeld L, Riecke L, Valente G, et al. Cortical Tracking of Multiple Streams outside the Focus of Attention in Naturalistic Auditory scenes. NeuroImage. 2018;181:617–626.
16. Puvvada, KC, Simon JZ. Cortical Representations of Speech in a Multi-talker Auditory Scene. J Neurosci. 2017;37(38):9189–9196.
17. Beck DL. Cognition, Audition and Amplification. Keynote presentation to the Missouri Hearing Society. Saint Louis, Mo. 2020.
18. Miller SE, Zhang Y. Neural Coding of Syllable-Final Fricatives with and without Hearing Aid Amplification. J Am Acad Audiol. 2021;31(8): 566-577.
19. Santurette S, Behrens T. The Audiology of More. Whitepaper. 2020.
20. Lesica NA. Why Do Hearing Aids Fail to Restore Normal Auditory Perception? Trends Neurosci. 2018;41(4):174–185.
21. Ronnberg J, Lunner T, Zekveld A, et. Al. The Ease of Language Understanding (ELU) Model: Theoretical Empirical and Clinical Advances. Front Syst Neurosci. 2013;7:31.
22. Santurette S, Ng EHN, Jensen JJ, et al. Oticon More Clinical Evidence. Whitepaper. 2020.