Plausible Northern Extension of Los Humeros Caldera Geothermal Field and Correlation of Density and Electrical Resistivity within the Collapsed Structure ()
1. Introduction
Los Humeros caldera is located on the eastern portion of the Trans Mexican Volcanic Belt (TMVB; Figure 1). This is a Pleistocene caldera with at least two eruption episodes that produced two nested calderas: Los Humeros and Los Potreros calderas; the latter is the youngest. Volcanic activity in Los Humeros caldera has been determined from 0.46 Ma to the Present [1] [2]. Large quantities of ignimbrites were produced in these occurrences. The whole area is underlain by limestone formations [3].
Figure 1. Digital elevation model (DEM) of the Los Humeros caldera (LH) region. Yellow triangles are monogenetic cones reported in the text. Major volcanic structures are: CP, Cofre de Perote; PO, Pico de Orizaba; LM, La Malinche. Contour interval 500 m from 1500 to 3500 m. Map from GeomapApp [4].
Geothermal studies at Los Humeros are directed at locating geothermal sources and those aiming to understand the structural characteristics of the caldera. The present generation of electrical energy at Los Humeros geothermal field exceeds 90 Mwe.
Los Humeros caldera is a superhot geothermal field where deep fluids reach temperatures of over 350˚C [5], location of the deep heat sources and its relation to the surface structures associated with caldera collapses are of vital importance for understanding the evolution of the geothermal system [6]-[9]. The aim of this study is to provide a new, high-resolution, gravimetric base map data of Los Humeros caldera, to better discern its internal structure to 6 km depth. The area covered with the high-resolution gravity is larger and denser than any of the previously reported gravity surveys. We explore the coincidences and discrepancies with geologic results, and two electrical resistivity studies obtained in magneto-telluric surveys performed in the caldera. Exploitation of this geothermal field has been executed for over 40 years with the drilling of numerous production and injection wells.
We use the gravity field of Los Humeros area obtained from the GGMplus gravity model [10] to create density models to depths of up to 6 km by means of 3D inversions, and extract density cross-sections from the volumes obtained. We compare our models with the geological section [3], and the electrical resistivity models obtained with the magneto-telluric method (MT) [11] [12], to find correlations between them.
1.1. Geological Studies
The first geologic descriptions of the Los Humeros caldera area appear to be those of Ordóñez [13] [14]. In the 1970s, the interest in the geothermal possibilities of Los Humeros caldera increased, and additional geological descriptions were made [15]-[17]; subsequent geologic and stratigraphic descriptions refined some of the concepts about the caldera [3] [18]-[21]. Isotopic and geochemical data were also included [22]. Accounts of the variations induced by exploitation of the geothermal field of Los Humeros are available [23].
1.2. Geophysical Studies
Reports on geophysical studies on Los Humeros caldera appeared almost 50 years later than the geological accounts; they start with a group of reports in a special issue of Geofísica Internacional [24] on Los Humeros caldera, covering studies on seismicity [25], gravity [26], magnetics [27], and tellurics, self-potential, and surface temperatures [28]. There are contemporary reports of Comisión Federal de Electricidad on electrical soundings, to locate potential geothermal sources [29].
More contemporary geophysical studies are available [11] [12], whose results we analyze here, and those that perform a joint inversion of gravity data and surface wave dispersion [30]. Their resolution is coarser than the resolutions we use here, since their inversions reach 20 km depth.
2. Satellite-Derived Gravity GGMplus
For the gravity field, we use the GGMplus gravity model [10]. The calculation of gravimetry in the GGMplus model is carried out by combining satellite, terrestrial, and topographic gravity data; the most recent gravity data measured by the GOCE and GRACE satellites are merged, offering information at scales of ~10000 to ~100 km. Data from the EGM2008 model is added, which covers scales from ~100 to ~10 km in the geoid [10] [31]. Finally, terrestrial, airborne, and marine data from the EGM2008 model are incorporated, combining them with gravimetry calculated from high-resolution topographic models, using the topographic gravity procedure, which is derived from high-resolution terrain models and complement scales from ~10 km to ~250 m [10]. This allows for the generation of the ultra-high-resolution gravity model with a spatial precision of 7.2 arc-seconds, and an accuracy ranging from 2 to 5 mGal [10], although this depends on the quality of the topographic model in the area and the gravimetric data collected for the model. In the Mexican highlands, we have determined an accuracy not exceeding 3 % error with respect to available land surveys.
![]()
Figure 2. Complete Bouguer anomaly map of Los Humeros caldera area. UTM coordinates correspond to Zone 14N. Gravity along lines A, B, and C, is extracted from the 3D gravity inversion. MT lines P6, P1, and P5 are obtained [11]. Line C corresponds to Line AA’ of Carrasco-Núñez et al. [3]. The dashed circle encloses the region of Los Humeros caldera. The regional gradient of the BA shows important alterations in the Los Humeros caldera region.
The Bouguer anomaly (BA) values in the BA map of Figure 2 show increasing values in the NE direction; those corresponding to Los Humeros caldera introduce a distortion in the gradient pattern; the area corresponds to the lowest gravity values of the BA map. To verify the relation of the surface gravity field with its potential sources, we performed a vertical derivative of the BA map. The result is shown in Figure 3, where the low on the northern rim of the caldera, and the smaller negative anomaly 12 km north are clearly identified.
2.1. 3D Inversion of the BA
3D gravity inversions were performed according to an established method [32], based on recognized theoretical considerations [33]; the Oasis Montaj program of Geosoft contains the inversion code. To represent geologic volumes, the program uses a Cartesian Cut Cell algorithm (CCC) to match the observed result with the calculated one; error limits can be expressed in terms of mGal, or in terms of a percentage of the standard deviation of the data (e.g., 5 %). The inversion program uses an Iterative Reweighting Inversion algorithm [34]. The inversion results are densities in g/cm3; in Figure 5, we show the observed and calculated gravity values along Lines A and B with a tolerance of 5% of the standard deviation of the data. The depth of the inverted volume is proportional to the size of the basic volume selected for the calculation; typically, choosing cubes 500 m on the side will reach depths of 5 - 6 km, and larger depths can be obtained at the expense of lesser resolution. After performing the inversion, one can explore the geologic volume obtained by extracting cross-sections, or density iso-surfaces that help visualize the associated structures.
![]()
Figure 3. First vertical derivative (Dz) of the Bouguer anomaly map in Figure 2 showing an N-S distribution of low values within Los Humeros caldera (dashed circle) with a minimum on its northern rim, which continues north revealing another, smaller minimum, not previously reported.
The potential discrepancies between the forward-modeled gravity data using an existing geological model and the observed gravity data, even after applying Butterworth filters, have been discussed [5]. Often, the geological model may not be able to represent the complex subsurface structures and variations in gravity anomalies. The authors conclude that the combination of high reactivation potential (i.e., resurgence) and low gravity along regionally significant structures that are NE-SW oriented in the study area indicates fracture porosity as the cause of the low gravity.
In the density contrast model [30], the same fault zones are characterized by comparably low density. To the east of Los Humeros geothermal field, density contrast is characterized by positive values. Los Potreros faults are associated with a density contrast of −400 kg∙m–3; such density contrasts along fault zones have been attributed to fracture porosity in other geothermal fields [35]. Here, we share the origin of the low-density areas between fractured porosity on the upper layers, and low-density magmatic materials in the lower layers.
2.2. Vertical Derivative of the BA
The vertical derivative (BA_Dz) of the Bouger anomaly is shown in Figure 4. This derivative enhances the location of the sources of maxima and minima of the BA map to which the derivative is applied. An N-S oriented, low-gravity anomaly is in the northern region of Los Humeros caldera, shown here crossing its northern limit. About 16 km north of the main anomaly, there is another well-delineated anomaly of shorter dimensions to be discussed ahead in relation to the density cross-section of Line A (Figure 2).
3. Density and Shallow MT Results
In this section, we present results from the 3D gravity inversion and complement them with results of a shallow MT study [12] restricted to the extent of LHC, that associates high- and low-density anomalies with results derived from the electrical distributions obtained with the MT study.
3.1. Line A
The density cross-section along Line A (Figure 4) shows that the three density regions approaching the surface are connected at depth, suggesting that they are fed by a single magma chamber located deeper than 6 km bsl. The association of low density and magmatic materials has been shown in several previous papers on volcanism: The intersections with Line B and MT lines P-6, P-1, and P-5 are indicated. The region corresponding to Los Humeros caldera is the region under the three MT lines; to the north there is an important low-density region that, up to now, has not been related to the caldera. Line B crosses this anomaly that we call the North Linear Extension of Los Humeros caldera, since it is the continuation of the N-S low-density regions within the caldera.
The density cross-section along Line A (Figure 2) extends 45 km in the N-S direction (Figure 4); intersections with the MT lines and Line B are indicated. The relevance of Line A is that it goes through three main gravity minima, two of them belong to LH caldera and the third one belongs to a region outside the caldera that appears to be associated with minor volcanic manifestations (Figure 1). On this extension, at (19˚49.11', 97˚25.60') we locate a monogenetic cone, near La Finca and Meyuco, and an additional one at (97˚27.20', 19˚48.16') SE from Tezhuatepec. Two additional cones may be located at (97˚24.46', 19˚49.62') and (97˚22.67', 19˚54.20'). Given the association with Los Humeros caldera, this region should be further explored for geothermal resources.
![]()
Figure 4. Density cross-section along Line A (Figure 2). Intersections with MT lines P-6, P-1, and P-5 [11], and Lines B and C are indicated. The density cross-section extends 45 km in the N-S direction and 9 km in the Z direction. Elements of the MT interpretation [12] are superposed to the density cross-section; their results extend vertically from sea level to +3000 m elevation. Isotherms, measured and inferred, are shown as red lines, from 150˚C to 250˚C in 50˚C intervals. LHC, los Humeros caldera. LHS, Los Humeros scarp. LPC, Los Potreros caldera. The scale represents density contrast values relative to +2.67 g/cm3.
Intersections along this density line with the MT lines indicate that P-5 and P-6 are located over negative density anomalies, whereas the central line P-1 is located over a positive density anomaly. Ahead, we will jointly analyze the density and the electrical resistivity responses to depths of 6 km bsl. Magma chamber location has been placed at depths of 5 - 9 km [6] [8]. The existence of a magma chamber extending from −3 to −9 km in depth, with two narrow regions at its edges, from the surface of the magma chamber to the base of LHC was proposed [36]; the authors analyzed thermal distributions with magma chambers at depths of 4, 5, and 6 km. In Figure 4, our density model confirms the existence of the two narrow low-density flow regions with a high-density region between them. Notice, however, that in the density distribution shown in Figure 4, the main chamber is not directly below the LHC. The main low-density region is located on the north half of the cross-section, exhibiting three density lows at elevations of +1000, +500, and −6000 m; we identify the latter as the magma chamber that feeds the upper deposits. The northernmost anomaly is the one proposed as a potential new geothermal area.
The central region shows the highest elevations of the caldera and a central high-density region surrounded by low-density conduits. We propose that the high-density region corresponds to the collapse that created the caldera, blocking the ascent of magmatic materials to the surface and inducing the bifurcation of the low-density regions around it.
Superposed to the density cross-section are the results of the MT survey [12]; they confine their results to the length of LHC horizontally, and vertically they restrict them from sea level to the surface, at +3000 m. In this projection, the elevated portion of LHC shows a central, high-density region flanked by low density anomalies that reach 9 km down from the surface. LPC is located on the S portion of LHC, between a scarp fault to the S and a normal fault to the N, whereas LHC continues up to the scarp fault to the N. Two arrows indicate a region of resurgent activity, mostly associated with the high-density region, in turn associated with a region of alterations caused by hydrothermal fluids containing iron and magnesium (propylitic alteration).
We reported on a similar density distribution in Tofua caldera [37], where the collapsed region of the caldera is of high-density, surrounded by low-density sectors; from that observation, we submit that a similar action may have occurred in LHC, where the collapsed portion is now experiencing resurgent activity in LPC.
3.2. Line B
Line B (Figure 6) intersects Line A perpendicularly at the location where the low-density anomaly appears more robust. Along this line, the width of the lowest-density region ranges from 3 to 5 km, indicating that injection of low-density material may be controlled by a normal fault extending N-S at least 30 km.
3.3. Line C
The density cross-section along Line C (Figure 2) was obtained exactly along the extent of Line AA’ of Carrasco-Núñez et al. [3] to compare their geologic interpretation with the density cross-section obtained from our 3D inversion. In Figure 7 we superposed the layers of materials (Qig, Tpa and K) reported by those authors to the density inversion results; we find that the thickest deposits of caldera components occur within a high-density region extending down to the –6 km elevation, or the bottom of the density model. Line A is located on the eastern rim of the high-density region, consequently, in Figure 4, it only samples a fraction of its volume. Line C shows a much larger version of this anomaly. It is to be expected that the deposits are deeper where the collapse was larger.
We submit that the configuration of the high-density anomaly to the W corresponds to an abandoned magma chamber, whose main deposit was located between sea level and –3 km elevation, as shown by the extent of the highest density region. After solidification, the weight of this region may have induced the collapse that created the caldera, coinciding now with the thickest deposits (Qig,
Figure 5. Observed and calculated gravity values along Lines A and B. The latter is the gravity field generated by the density distribution obtained from the 3D inversion. The error limits are established by the operator prior to the inversion process, either by a percentage of the standard deviation of the data (e.g., 5%, which is the case of the data shown), or by defined gravity values (e.g., 2 mGal). The greater the precision, the longer the computer time.
Figure 6. Density cross-section along Line B (Figure 2). The intersection with Line A is indicated. The region of lowest density occurs in a region ~4 km wide, and the full, low-density vertical region is ⁓8 km. The scale represents density contrast values relative to +2.67 g/cm3.
Tpa, K) within the caldera. We have also observed that the high-density regions displace the original ascent trajectory of the low-density ones, which appears to be the case of the central low-density anomaly. We assume that volcanic activity initiated along with the present high-density region, which eventually became obstructed, stopping magmatic supply along that trajectory; cooling of that region increased its density, eventually inducing the collapse, which originated the caldera structure.
Magmatic fluids could not flow through the high-density region, being displaced to the east, forming the central, low-density conduit.
We refer to the low-density distribution as a feeding channel since we have repeatedly observed similar low-density distributions associated with active volcanoes (e.g., Colima Volcanic Complex [38] [39]; Popocatépetl [40], and volcanic calderas [37] [41]). Furthermore, the easternmost low-density region is interpreted as a conduct, like the central one, but with minor, or no surface manifestations. Notice that both low-density channels present lows at −6 km, which appear to be associated with magmatic chambers fed by deeper magmatic sources.
![]()
Figure 7. Density cross-section over Line C, Line AA’ of Carrasco-Núñez et al. [3] (Figure 2). The interpretation of these authors, of the internal structure of the caldera from +1000 m to the surface, is superposed on the density cross-section. The depth distribution of layers of limestone rocks, andesites and dacitic lavas, and basaltic materials are shown (magenta). Elements of the MT interpretation [12] are superposed to the density cross-section; their results extend vertically from sea level to +3000 m elevation and horizontally extend only the length of Los Humeros caldera. Isotherms of 150˚C, 200˚C, and 250˚C are shown in yellow. LHC, Los Humeros caldera. LPC, Los Potreros caldera. OS, Oyameles scarp. LHS, Los Humeros scarp. The scale represents density contrast values relative to +2.67 g/cm3.
Propylitic alterations occur more often here than in the N-S projection (Figure 4), being larger in association with the high-density region, possibly indicating, in conjunction with results of the N-S projection, larger hydrothermal fluid fluxes in the more fractured, high-density resurgent region. Also, in this projection, LPC corresponds with the high-density region. In this projection, the high-density region extends westward of the Oyameles scarp, exhibiting a maximum density value centered at an elevation of −1000 m, directly connected to the W limit of LH and LP calderas, also showing various wells drilled in the region, and overlapping the largest propylitic alteration area. The E limit of LHC corresponds to a low-density channel with a propylitic alteration close to the surface, indicating this region also contains relevant hydrothermal fluid flow. The MT study [12] does not reach the easternmost, low-density anomaly; however, given its similarity with that of the eastern rim of LHC, some geothermal exploration in the area appears advisable.
We note the similarity between cross-sections B and C; considering their separation of over 15 km, we find that the N-S structural characteristics extend beyond the north rim of the caldera, particularly the presence of the thin, low-density regions in both. From the cross-sections along Lines A, B, and C, we infer that magmatic material has been fed through a long (~30 km), ~13 km wide section, containing a narrow (~4 km) conduit of lowest density that extends north of the caldera limits. This interpretation is validated by the result of the vertical derivative (Figure 4) which confirms the extent of the low-density anomaly north of the caldera rim.
4. Density and Electrical Resistivity Comparisons
Magneto-telluric (MT) soundings were performed in Los Humeros caldera along Lines P-6, P-1, and P-5 (Figure 2) [11]. We extracted density cross-sections along those lines at a model resolution of 500 m (i.e., prisms with side dimensions of 500 m); this resolution was chosen to reach depths of 6 km bsl, to allow an adequate comparison with the depths obtained with the MT method.
4.1. Cross-Section P-6
Figure 8 shows the density cross-section with superposed results from the MT P-6 resistivity cross-section [11], the location of two geothermal wells, and ten MT stations. We first describe the characteristics of the density cross-section. A low-density, vertical region appears between X = −5000 and +5000 m, the central portion contains the lowest density values and on the extremes two low-density regions rise toward the surface. The lowest values of the central, low-density anomaly are located between −5000 and +1000 m elevation. We consider this plane intersects the magmatic deposit between −5000 and −1000 m, where the lowest density magmatic products are located. The small region located under MT stations 29-32 is interpreted as a hot intrusion of low-density materials close to the surface, which must be contributing to the heating of that surficial volume; it is in this region defined by the brown, horizontal layers in Figure 8, that low resistivities are located.
Superposing the electrical resistivity results, we observe that between the two HR regions, there is a vertical channel (IR) where resistivity decreases. We note that it corresponds to the region of lowest density, leading us to infer that this more electrically conductive path corresponds with the presence of melt in the
![]()
Figure 8. The density cross-section corresponding to the trace of the MT profile P-6 [11] extends along 30 km on Los Humeros caldera (Figure 2); the vertical yellow lines delimit its extent on the larger density cross-section. The numbers in blue correspond to MT stations. HR denotes high resistivity (>1000 ohm-m) and IR denotes intermediate resistivity values (150 - 300 ohm-m) constrained between yellow lines. Low resistivity (LR) is <50 ohm-m. In addition to the surface layer (black), two brown layers define the region where low resistivities are located with the MT soundings. Faults are indicated by dashed lines and the depths reached by the wells correspond to straight lines. UTM and geographic coordinates appear above. The scale represents density contrast values relative to +2.67 g/cm3.
magmatic chamber; it involves an area of ⁓12 km2. HR denotes high resistivity (>1000 ohm-m), IR denotes intermediate resistivity values (150 - 300 ohm-m), and LR refers to low-resistivity (<50 ohm-m). High resistivities tend to follow the two lateral low-density regions that rise toward the surface. We note that the low-density region may be considered as a vertical channel ⁓10 km wide from the surface (+3000 m) to the bottom of the cross-section (−6000 m) that is properly represented by the high-resistivity regions. The top of both HR regions adequately conforms to the low-density distributions above sea level. The averaged elevations of the gravity inversion appear as a raised, central section that corresponds to the region of low densities below. The raised portion corresponds to higher densities than those of the substrate, and lower resistivities; the former may arise from the cooled components of the caldera collapse, while the low resistivities may originate in the presence of aquifers, from which the geothermal resource originates.
As in the case of Line A (Figure 4), the central, high-density region is flanked by low- density regions; this configuration extends down at least to −1 km [42] registered simultaneous eruptions at the E-W caldera limits at 7.3 ka at the Cuicuiltic member eruption, which can be directly connected with these low-density regions, indicating an active role of the low-density regions in the volcanic processes involved in LHC and LPC.
Summarizing the results of the density-electrical resistivity correlation along profile P-6, we first note the low-density, vertical region between X = −5000 and +5000 m dominating the central portion of the cross-section. Next, we underline the correlation of a central channel of intermediate resistivity that overlaps the lowest density region, and finally, we also note the coincidence of the surficial intrusion of mid-density material and the low-resistivity regions found with the MT method.
4.2. Cross-Section P-1
Cross-section P-1 runs E-W at the middle of the caldera (Figure 2). The cross-section is divided into a high-density region to the W, and a low-density section to the E (Figure 9); both reach the surface. The corresponding portion of Line C (Figure 7) shows a similar density distribution in the portion that can be compared; the density distribution directly along P1was obtained and cannot be distinguished from that obtained for Line C. The electrical resistivity distribution shows regions of high resistivity (HR), intermediate resistivity (IR), and low resistivity (LR). A channel of low to intermediate resistivity is present on the western portion of the section, from depths of −6 to +2 km, which ends on the surface in a small region of intermediate density, distinct from the higher values surrounding it. Contrasting with the observation of a similar channel in cross-section P-6, this channel occurs in a high-density region. We speculate that the high-density region corresponds to the resurgent portion, or LPC, where high porosity has developed, allowing for a larger geothermal flux and consequently lowering the electrical resistivity along that channel. In this case, uplift of the high-density block would be produced by new, low-density magma pushing the block upwards.
The purple region starting on the surface (Figure 9) and deepening towards the E corresponds to a group of low-resistivity sections; it correlates with regions where the high-density decreases. The LR region close to the surface deepens towards the E, crossing the sharp, high-low density boundary. Notice that Line C and P-1 are not along the same orientation; there is an 11˚ angle between them (Figure 2). After obtaining the cross-section corresponding exactly to P-1, we observed that both cross-sections are quite similar, maintaining the original Figure 9.
We interpret the upper LR region as one where low-density materials have been reaching the surface, probably since 0.46 Ma [1] [2]; the high-density region may have been intruded first as a low-density region, magma feeding at this region stopped, starting to cool off and becoming a high-density one, while the low-density one remained, maintaining hotter, thus less dense, magmatic materials. Since the high-density region corresponds to the resurgent area, the associated motions may have induced a larger porosity. Although speculative, these statements offer a potential sequence of events in the region.
![]()
Figure 9. The density cross-section corresponding to the trace of the MT profile P-1 [11] extends along 25 km on the central portion of Los Humeros caldera (Figure 2). There are regions of high electrical resistivity (HR), intermediate resistivity (IR), and low resistivity (LR) whose values in ohm-m appear in Figure 8. The numbers in blue correspond to MT stations. There are four geothermal wells. The violet region close to the surface corresponds to a series of low-resistivity sections that are lumped together. The scale represents density contrast values relative to +2.67 g/cm3.
4.3. Cross-Section P-5
The northernmost cross-section (Figure 2) runs parallel to P-1 and P-6, with its eastern part extending in a smaller line at a different orientation. The density cross-section P-5 does not include the latter extension of the line. Figure 10 shows the superposition of the electrical resistivity regions on the density cross-section.
Figure 10 shows the density cross-section with superposed results from the MT profile P-5 resistivity cross-section [11]. The intermediate resistivity channel observed in P-6 and P-1 appears again in this projection; however, its location is close but displaced to the west from the region of minimum density.
4.4. Resolution Comparison
We present the comparison between inversion resolutions of 500 and 250 m (Figure 11) to illustrate the importance of performing inversions at higher resolutions. In principle, the higher the resolution, the closer the model to the actual density distribution; the density values of the model also change accordingly. Although
![]()
Figure 10. Regions obtained from the MT profile P-5 (Figure 2) are superposed on the corresponding density cross section. HR denotes high resistivity (>1000 ohm-m), IR denotes intermediate resistivity values (150 - 300 ohm-m) constrained between yellow lines, and LR denotes low resistivity (<50 ohm-m). In addition to the surface layer (black), a brown layer defines the region where low resistivities are positioned with the MT soundings. Faults are indicated by dashed lines and the depths reached by the wells correspond to straight lines. Numbers in blue correspond to MT stations. UTM and geographic coordinates appear above. The color scale represents density contrast values relative to +2.67 g/cm3.
![]()
Figure 11. Two versions of the density cross-section P-5 are presented, the one with a resolution of 500 m is the same as that in Figure 10; superposed to it is the one at a resolution of 250 m. Their horizontal scales are matched, while the vertical scale of the latter is displaced downwards to allow for the comparison of results between +3000 and -1000 m. The scale represents density contrast values relative to +2.67 g/cm3. The scale of the superposed cross-section ranges from −0.20 to +0.07 g/cm3.
the local density values change, the trends are maintained. In Figure 11, we can see how the large negative anomaly at 500-m resolution, under MT stations 22 - 27 becomes a fragmented, negative channel in the 250-m resolution version, where only the uppermost sections contain the lowest-density values. Similarly occurs with the high-density region under MT stations 19 - 21. In one case, the volumes’ density averages involve cubes of 500 m on the side, whilst in the second case, the averages are on cubes of 250 m on the side, which allows for a better representation of the local, ground variations. Modeling larger depths is possible at the cost of lesser resolution; the case of the Popocatépetl volcano is at hand, where model depths of 17 km were attained with a resolution of 2.5 km, to compare with geothermometry results [40]. We note that in the case of MT surveys, the resolution of the measurements cannot be so easily modified, since resolution depends on the sounding frequency, and the station separation. We used profile P-5 to illustrate the resolution effects on the inverted data, but similar results at the 250-m resolution were obtained for profiles P-1 and P-6, which are not included for reasons of space. An illustrative example of the resolution role on the representation of an anomaly can be found in a recently published article on the magmatic chambers of Izaccíhuatl volcano [40].
5. Discussion
The perpendicular density distributions (Figure 2) in Figure 5 (N-S) and Figure 6, Figure 7 (E-W) indicate that the magmatic deposit under LH caldera is elongated in the N-S direction but short in the E-W direction, suggesting that the deposit has been placed using a normal, geologic fault running in the N-S direction. The density distribution in Figure 7 that corresponds to Line AA’ of Carrasco-Núñez et al. [3] shows a good correspondence between their layered materials’ distribution and the density distribution obtained from our 3D inversion model. We infer that the high-density region surfacing between 660000 - 665000 N (Figure 7), corresponds to an older intrusion, today forming part of the resurgent section. Regarding the statement that the magma chamber presently feeding Los Humeros caldera is displaced 12 km from the caldera center, we obtain the first indication from the vertical derivative (Figure 3), where the main minimum, elongated in the N-S direction, is in that position. This evidence is confirmed in the N-S density cross-section (Figure 4), where the density minimum at 9 km depth (6 km bsl) appears in that position feeding Los Humeros caldera and a second branch surfacing at 19˚50'N, which is the newly proposed region with geothermal potential.
We observed a vertical, intermediate resistivity channel in the three electrical resistivity sections analyzed. The intention was to associate it in every instance with the low-density regions; however, we think that the association is better made with regions in which bulk porosity favors geothermal circulation inducing propylitic alterations.
We find good agreement with the results [30] along their Line AA’ and our Line A (N-S), and their Line BB’ and our Line B (E-W). Their lines reach 20 km depth, while ours reach only 6 km, showing more detailed density distributions.
We have interpreted the central high-density region as induced by the collapse of the central portion of the LHC, mostly based on observations made in other caldera-forming, volcanic structures [37] [40]. The association of low-density regions with magmatic materials in volcanic settings is well established [5] [30] [37]-[40]. Often, the vertical distribution of these features is an additional indication of their relationship with magmatic materials. Under these circumstances, variations in basement lithology or the presence of thick, low-density sedimentary units tend to be considered a minor possibility.
6. Conclusions
Location of the magmatic chamber responsible for the geothermal activity in a given area is of primary importance; since this parameter is fixed, modeling of its geological or geophysical characteristics results in more constrained exercises. Contrasting with most of the proposed locations in the literature, the magmatic chamber of LHC is not directly located under LHC; our results indicate it is displaced ⁓7 km to the N, with one arm of the low-density feeding trajectories located on the S portion of LHC.
Having established that the magma chamber is displaced N of the location of LHC, we located a new, low-density region connected to the same magma chamber that has not been reported previously. This region reaches the surface, as does the LHC low-density region. Given its interconnection with the geothermal field of LHC, as well as the presence of some smaller volcanic structures, we suggest carrying out additional geothermal studies in this area, such as drilling exploratory wells, geochemical analyses, and magneto-telluric and seismic surveys. Detection of the continuation to the north of the low-density anomaly of Los Humeros caldera is a significant addition to the regional characterization of this important geothermal region.
Based on geological and geophysical evidence, including the presence of calcite-wollanstonite in wells H-6, H-8, H-9 and H-12, probably caused by a metamorphic contact with Cretaceous limestones, a cylindrical magma chamber with thin vertical conductors at its extremes was proposed [36]. Our results along cross-sections of Line A, P-5, and P-6 show the two vertical, low-density regions, which coincide with the peculiar characteristics of the magma chamber proposed by the above authors.
From the electrical resistivity comparison with density distribution, we conclude that IR and LR occur associated with either highly fractured regions with larger porosities, or low-density regions where magmatic materials accumulate. The low-resistivity regions were found on the upper portion of the MT sections [11], coinciding with the uppermost layers of the caldera formations, which in turn correspond to regions in which rapid variations of density were observed. These observations are not enough to pinpoint specific locations for the exploitation of the geothermal resource; however, they contribute to the general understanding of the complex characteristics of the geothermal reserve in Los Humeros caldera.
Acknowledgements
This study has been supported by IIMAS and Instituto de Geofísica, UNAM; we acknowledge material support from both institutions. This research did not receive any specific grant from funding agencies in the public, commercial, or not-for-profit sectors.