Southward Migration of Magmatic Activity in the Colima Volcanic Complex, Mexico: An Ongoing Process

The Colima Volcanic Complex trends in a nearly N-S direction in western Mexico, and one of its structures, Colima volcano, is the most historically active volcano in the country. Immediately to the N, there is another volcanic center called El Cantaro volcano, whose activity started around 1.7 Ma in its N portion and migrated to the S in various episodes. Volcanic activity migrated further south, from El Cantaro to the Colima Volcanic Complex where the southernmost manifestation, Hijos del Volcan domes, is located on the south slope of Fuego volcano. The above date appears to mark initiation of the rather continuous volcanic activity in the area. It has been noted that these volcanic manifestations lie on, or near the Rivera-Cocos inland plate boundary. Colima’s Fuego volcano is also the closest to the Middle America Trench, among the polygenetic volcanoes in Mexico. We submit that the anomalous location of volcanism in this area originates in an anomalous subduction process of the Rivera and Cocos plates and evoke a tectonic model, proposed elsewhere, to support the idea. Modeling gravimetric and aeromagnetic data we locate the magma chambers of the Fuego (active) and Nevado (extinct) volcanoes within a 65 mGals negative Bouguer anomaly elongated in a nearly N-S direction. The corresponding aeromagnetic map displays a magnetic high over the southern portion of the Fuego volcano edifice. We found two additional, associated structures whose anomalies have not been previously reported, which appear to follow the southward magmatic migration pattern. One of them is a collapse structure with a circular topographic expression, and the southernmost is a low-density intrusion ~1 km below sea level, associated with a moderate topographic bulge at the surface that we interpret as a magma body. Five lines cross the anomalies; gravimetric and magnetic fields are concurrently modeled along them to locate the magmatic bodies. In addition to the 2-D models we perform 3-D gravimetric and magnetic inversions. For each field a 3-D mesh is built under the area occupied by the Colima Volcanic Complex, the volume elements are then assigned density or magnetic susceptibility values and their surface contributions in various points are evaluated. The process is iterated until the difference between the measured and the calculated fields is less than a predetermined value. The results of each inversion adequately and independently define the location of the magmatic chambers although they cannot distinguish between the individual chambers of the Nevado and Fuego volcanoes. 2-D and 3-D results complement each other and consistently show the locations of potential magmatic regions. Our models support a multiple, complex magmatic system that appears to continue to spread southwardly, which can pose additional volcanic risks to an already threatened local population.


Introduction
The term Colima Volcanic Complex (CVC) was coined by [1] and comprises the active Colima Fuego (3860 m) volcano as well as the extinct Colima Nevado (4340 m) volcano. These structures are located in the vicinity of (19.54˚N, 103.61˚W) in western Mexico. To the N of this point there is, however, an older volcanic structure: El Cántaro volcano, with an N-S elongated, dissected structure around (19.72˚N, 103.63˚W). Figure 1 shows the location of these volcanic edifices as well as the general tectonic features framing them. In [1] they also noted that the southward migrating volcanism from Cántaro to Volcán de Fuego chain is also reflected by a southward-younging progression at Cántaro itself, giving rise to the two volcanic edifices of Cántaro volcano. Reference [2] placed the Potassium-Argon ages of Cántaro from 1.7 to 1.0 Ma. Volcán de Fuego is presently the most active volcano in Mexico, thus framing the rather continuous volcanic activity in this area between 1.7 Ma and the Present. Volcán de Fuego de Colima is often simply referred to as Colima volcano and the extinct one as Volcán Nevado de Colima or Nevado volcano.
The CVC is anomalously located with respect to the Trans-Mexican Volcanic Belt occupying the central portion of the Colima rift. Reference [3] described the Colima rift, or graben, as a N-S-oriented structural depression that coincides with a 100-km offset in the volcanic front, making the Colima volcano ( Figure 1) the closest to the Middle America Trench (150 km). Its existence has been connected with a transtensive regime originating in the separation of the Rivera and Cocos plates at depth [4]- [6]; however, the mechanisms inducing the volcanic activity and its southward migration have not been explained. The question arises whether such migration has stopped or is still occurring and where.
In this work we present gravimetric and aeromagnetic data in the region of the CVC, as well as the associated anomaly maps, several 2-D models along selected cross-sections, and two inversions, one gravimetric and the other magnetic, of such fields. The CVC is characterized by a broad negative gravimetric anomaly that encompasses the two volcanos: the active (Fuego) and the extinct (Nevado) structures. Negative volcano anomalies are usually associated with silicic volcanism [7]. Reference [8] presented a model of the magma chamber associated with the Colima volcano based only on gravimetric data, while [9] presented a model that includes Nevado and Fuego volcanoes based only on aeromagnetic data; the latter author emphasized the difficulty of modeling the aeromagnetic data since the magmatic body is assumed to be above the Curie point and thus, devoid of magnetization. It is interesting noting that the latter author's model [9] shows a thin (−6500 < z < −5500 m) extension of the magma chamber of the Fuego volcano, reaching ~6.5 km towards the south of its chimney. This result is refined in the 2-D models that we present here, which include the combination of both gravimetric and aeromagnetic data and are thus better constrained than the previous models. The baseline along the CVC is not only larger than those in the previous models, but is also crossed at selected points with transverse lines that complement and further restrict the present models. The gravimetric and magnetic field inversions provide 3-D views of the underground density and magnetic susceptibility distributions that help visualize the extent and location of  [34] triple junction (RRR, yellow symbol). To the S of the CVC there is a region often designated as the Southern Colima Rift (SCR), where the rift structure has not, however, been identified. Continuing S and offshore, the Manzanillo Canyons are shown as thin lines (Manzanillo Canyon MC, Cuyutlán Canyon CC, Armería Canyon AC, Coahuayana Canyon CO). The dashed line starting at the trench and ending at the CVC is the Rivera-Cocos plate boundary proposed by [6]. In the offshore region it virtually superposes to the Armería Canyon. The Middle America Trench (MAT) is the dark blue region that, offshore, follows a course similar to that of the coastline. The Rivera, Cocos, and North America (NAM) plates, as well as locations of the Rivera Fault Zone (RFZ), Manzanillo (M), Puerto Vallarta (PV), Bahía de Banderas (BB), the Jalisco Block (JB) and its NW limit (NWL) as proposed by [35], are also identified. magmatic bodies.

Gravimetry
A 260 gravimetric station data set is shown in the volcanic area of the CVC (Figure 2). The map shown in this figure corresponds to the residual Bouguer anomaly, obtained from the total, complete Bouguer anomaly subtracting the regional trend with a second-degree polynomial. Stations were measured with a Scintrex Autograv (CG5) gravity meter with a reading precision of 1 µgal; latitude, tidal, drift, and topographic corrections were made. The complete Bouguer anomaly was calculated using a reduction density of 2.67 g/cm 3 ; the estimated total error is below ±0.30 mgals. In order to obtain a uniformly spaced mesh the set was interpolated with a minimum  curvature algorithm that generated a 400 m grid. Contours are shown for 0 to −35 mGals in 5 mGal intervals, although the complete anomaly displays a ~65 mGal amplitude. The position and elevation of each station were obtained by means of a Garmin GPS with estimated average errors of 5 m in the X-Y plane and 8 m in the vertical direction. The vertical uncertainty is the most critical; taking into account the gravity gradient of 0.3086 mGals/m, we have as a maximum error of gravity anomaly up to 2.46 mGals. Since the range of anomalies observed in the study area is about 20 -30 mGals, the uncertainty in the height determination would induce a maximum uncertainty in our estimates of Bouguer anomaly of 8.2 percent, which is acceptable in regional gravity studies. Station distribution was determined by road and trail accessibility; the influence of interpolation artifacts, if any, is avoided analyzing only gravimetric anomalies of wavelengths several fold the separation between lines or stations.
The prominent gravimetric low delimited by the −35 mGal contour actually consists of two neighboring lows divided by a WNW-ESE trending high; this high is well defined by the gravity station distribution in the area. The larger anomaly to the N corresponds to the Nevado and Fuego volcanoes, while to the S the gravimetric low, not previously reported, may correspond to another region of magmatic activity, as suggested by the models to be presented ahead. To the N of these two gravimetric lows lie the southward-younging structures of El Cántaro volcano [1]. We will refer to the two structures as Cántaro South (CS) and Cántaro North (CN) based on topographic evidence and the different ages reported for them (1.7 Ma for CN and 1.0 Ma for CS [2]). Although this study does not intend to model Cántaro structures it is worth mentioning that these structures lie on a gravimetric low of lesser magnitude than that corresponding to the Nevado and Fuego volcanoes and is aligned with their main gravimetric anomaly.

Aeromagnetics
The aeromagnetic map used is that in [10], which we reduced to the pole, a technique that shifts magnetic anomalies to a position directly above their sources [11], and upward continued it [12] to 2 km in order to smooth the short-wavelength variations present closer to the terrain surface. The grid generated was a 400 m uniformly spaced grid shown in Figure 3. The general trends are quite similar to those in the aeromagnetic map of [9], as could be expected.

2-D Models
We calculated the gravity and aeromagnetic models along five lines that sample the anomalies corresponding to a roughly SW-NE profile (L0, azimuth 28˚) and four perpendicular profiles (L1 to L4) to this one (Figure 2) crossing the main structures in the area. A general description of the Bouguer anomaly map (Figure 2) is in order at this stage. The anomaly's amplitude is 65 mGals with a wavelength of 50 km, fitting well in [7] classification of silicic volcanic structures. In the map the anomaly extends along L0 and is bounded by the −35 mGal contour. The gravimetric high dividing the anomaly, slightly south of L3, decreases its amplitude towards the SE and suggests the existence of a channel connecting the northern and southern portions of the split anomaly. This channel aligns with the Nevado, Fuego, and Hijos del Volcán structures. Regarding the aeromagnetic anomaly map in Figure 3 we can point out in a general way that it resembles a dipolar anomaly with a wide positive area to the south and a smaller, negative one to the north, coinciding fairly well with the gravimetric low bounded by the −35 mGal contour. The largest positive amplitude of the anomaly is located between Fuego volcano and Hijos del Volcán domes.

L1 Nevado Volcano
The model of the Nevado volcano along L1 is shown in Figure 4; it has a horizontal extent of 18 km and is modeled down to 9 km depth. The topography along this line and all others modeled here is obtained from the 30-m digital elevation model of [13]. The amplitude of the residual Bouguer anomaly along this line is 32 mGals and two magma chambers at different depths are represented by bodies of 2.35 g/cm 3 density value with negligible magnetization ( Table 1), suggesting that, although the volcano is extinct the magma chambers may still be above the Curie point. The aeromagnetic anomaly with approximate amplitude of 200 nT is shown in the uppermost portion of the figure together with the calculated model. The length of the upper magma chamber along this direction is 6 km and is located at a depth of ~1.2 km below sea level; its average thickness is ~3 km. The remains of the volcanic chimney merge with other volcanoclastic formations of densities 2.58 to 2.65 g/cm 3 above the upper magma chamber. Geological formations with densities between 2.75 to 2.80 g/cm 3 underlay the  [10] reduced to the pole and upward continued 2 km. The grid generated is a 400 m uniformly spaced grid. The gravimetric contours are superposed, as well as the modeled lines (L0-L4); abbreviations are the same as in Figure 1. A magnetic high is associated with Colima Fuego edifice extending southward to the vicinity of the proposed magma body (M). A magnetic low is located between Nevado and Cántaro Sur structures, coinciding with the north limit of the gravimetric minimum. Magnetic minima continue northward along the northern Colima graben. To the W another deep magnetic low is present. Coordinates are geographic and NAD 27/UTM Zone 13N. magma chamber. The model also suggests the existence of a smaller magma chamber located at a depth of ~5 kmbsl, directly below the upper magma chamber. The larger magma chamber appears to be enclosed in the Upper Cretaceous limestone that outcrops west and east of the CVC [14] while the smaller magma body is overlaid . E-W gravimetric and aeromagnetic model along L1, Nevado de Colima. Two, superposed magma chambers are shown; the top of the shallower one is at ~1.2 km below sea level, shown in pink to denote an assumed temperature lower than that of the deeper, smaller magma chamber. Both chambers are modeled with no magnetization, implying that both are above the Curie point temperature. The larger one would have fed the eruptions of the primitive Nevado volcano and would be in the process of cooling and crystallization. Several dashed lines labeled C1, C2, and C3 represent the intersection of L1 with three caldera rims mapped by [15]; C1 is the oldest and C3 the youngest. C2 does not have a counterpart in the east section. The green formation represents Upper Cretaceous limestone overlaid by volcanoclastic deposits.
by the limestone and underlain by the formation of density 2.8 g/cm 3 .
In [15] they reported at least three major collapses of this volcano during its active period. Reference [16] mapped three different calderas that resulted from the collapses. The locations of the rims of the collapsed calderas (C1, C2, and C3) and the location of the associated normal faults (dashed lines) are also shown in Figure  4. Arrows indicating the relative motions between the blocks are shown only in one line in order to avoid cluttering the figure. C2 does not have a counterpart in the east portion of the model probably owing to its destruction by the subsequent formation of C3. Only the eastern rim of C1 apparently lies outside the limits of the upper magma chamber, although it could be associated with the limit of the lower magma chamber.
The regional N-S trending fault system associated with the northern Colima rift is well documented [3] [17] [18]. In addition, [19] reported a new, transverse fault system across Nevado volcano. These two orthogonal fault systems intersect in the CVC area. The E-W trending fault system creates a shallow graben centered beneath the summit of Nevado de Colima. According to [20] graben structures occur in volcanoes upon magma intrusion. Their models show that the intrusive complex develops in continual competition between upward bulging and lateral gravity spreading. They define volcano spreading as the lateral deformation of the volcano edifice Table 1. Physical properties of the various layers used in the models. Mixed volcanic deposits contain airflow deposits and lava flows. Volcanic deposits contain lava flows and debris avalanches. Depth is positive downwards, as in the models. See  along a basal ductile layer that deforms under the volcano weight. Grabens, en echelon faults, folds, and thrusts are the characteristic associated structures. Varying parameters such as height/radius of volcano, or volcano height/substratum thickness they obtained the associated formation of single or multiple grabens. For brittle layer/ductile layer ratios >1.25 they obtained a single graben that traverses the volcanic edifice and that passes above the intrusive complex, resembling the case of Nevado volcano structure. Reference [20] validated their model results with field observations on eroded volcanic remnants as well as on active volcanoes by means of geophysical evidence. They offer the example of a dormant volcano, Piton des Neiges of La Réunion Island in the Indian Ocean, where a gravimetric survey shows a close association between the (positive) gravity anomaly and the star-shaped pattern of escarpments on the volcano. They propose that the graben structures and the associated fracturing developed are linked to the spreading of the edifice and to the growth of the intrusive complex. Similarly, we propose that the graben structure in Nevado volcano developed under the emplacement of the body described here as the upper magma chamber. Note, however, that Nevado's structure involves a negative gravity anomaly.

L2 Fuego Volcano
In Figure 5 we show the model of the magma chamber of the Fuego volcano (2.31 g/cm 3 density); this line goes exactly through the summit of the volcano. The magma chamber along this projection shows two lobes, both located between 1 km above and 2 km below sea level; the eastern lobe has a radius of ~3 km, matching the dimensions obtained for the chamber along the perpendicular line L0, to be presented ahead, while the western lobe is smaller, of approximately 1 km in diameter. This model differs from that of [9] since the corresponding vertical dimension in his model is almost 7 km. The magnetization of these magma bodies is also negligible since in this erupting edifice they must be above the Curie point temperature. A 2.58 g/cm 3 density conduit represents the volcanic chimney. The density of the geologic formations above the chamber (pyroclastics and andesitic lava flows) range from 2.55 to 2.6 g/cm 3 , while the Cretaceous limestone formations fully contain the magma chamber. The assigned density of the Cretaceous limestone is 2.75 g/cm 3 and that of the formation underneath them is 2.8 g/cm 3 . Two dashed lines show the intersections of the collapsed caldera of the ancestral Fuego volcano [21], a structure that produced a large debris avalanche towards the south dated by these authors at 4300 years B.P. The topography of the caldera rim is neatly shown on the western side while the eastern one has been obliterated by the lava flows of the new volcanic edifice that started growing around 4000 years ago [14] [21]. As in the case of Nevado, there is also a good correlation between the dimensions of the modeled magma chamber and those of the collapsed caldera of the ancestral Fuego volcano.

L3 La Escondida
Figure 6(a) shows a 3-D view of La Escondida and Fuego volcano where a depression is highlighted and some of the large hummocks inside it can be appreciated. In Figure 2 the location of the perimeter of this structure appears around the intersection of L0 and L3. Various authors have made reference to this structure. Reference The magma chamber eastern lobe radius is ∼3 km and the western lobe radius is ~1 km. The intersection of L0 with this line is shown, as well as two dashed lines that represent the intersections with the rim of the caldera of the ancestral Fuego volcano. There is a good correlation between the size of the collapsed caldera and that of the modeled magma chamber.
[21] described it as an area of ~20 km 2 , E of San Antonio Ranch, planned to a S-SW sloping surface where isolated hummocks project up to several hundred meters. They showed the area in a topographic map contoured at 100-m interval, but made no reference to its circular shape and no name was assigned to it. The depression was attributed to erosion of the area. Reference [16] described this region as composed of lahar-like deposits, interpreted as reworked hummocks of the debris-avalanche deposits, but made no reference to its circular shape or identified it as a caved in region. Reference [22] first described it as a collapse structure, observing its circular surface expression and pointing out the potentially high volcanic risk character of the structure. No explanation for the existence of the depressed or collapsed region was offered. These authors mapped the extent and geometry of the structure from satellite images and named it La Escondida (the Hidden One) after one of the lakes in its southern, internal margin. Reference [23] referred to this structure as El Jabalí depression after another lake neighboring La Escondida Lake. They associated it with a tectonic origin, but did not show a map of the structure. Reference [19] showed a map of the structure on a satellite image and called it La Yerbabuena depression after a small village in its internal, NW portion. They argue that this is not a tectonic feature, joining the opinion of [21] that attributes its origin to an intense erosion process. Considering the variety of names given to the structure we will follow the customary practice of adhering to the first published name: La Escondida collapse or depression. In Figure 6(b) we show La Escondida depression and three lines crossing it at three azimuths that intersect each other at the center of the structure. Figure 6(c) shows the topographic profiles along these lines, where arrows indicate the intersections with the collapsed perimeter of the structure. Vertical displacements from 50 to 200 m are observed around the perimeter, while in the southern rim there are several small lakes (La Escondida, La María, El Jabalí) often present in the lower topographic areas of collapsed structures. The general slope of the structure in the SW-NE direction is 5 degrees while in the NW-SE direction is 2 degrees. The age of the avalanche deposits is a matter of debate; [16] dated the pyroclastic layer overlying the avalanche deposits of the ancestral Colima volcano at 9370 ± 400 yr B.P. while [21] reported a 14 C age of 4280 ± 110 years. Regardless of this age discrepancy one can be assured that the collapse postdates the avalanche, since the interior of the collapsed structure shows avalanche deposits profusely; that is, the collapse occurred after the avalanche produced by the ancestral Colima volcano. In the surface a 7 km diameter circle can approximate the perimeter of the structure. Figure 7 shows the model corresponding to La Escondida structure. The modeled, tabular intrusion (sill) has a thickness of ~1.5 km and is located at a depth of ~1.2 km below sea level, coinciding with the magma chamber depths of neighboring Nevado and Fuego volcanoes. We submit that, before the structure collapsed, a magma Figure 7. E-W model along L3 La. Escondida structure. The 2.78 g/cm 3 density intrusive may have been injected molten as a sill. Subsequently magma injection stopped and the body started to cool down increasing its density and inducing the collapse observed in the surface whose perimeter appears in Figure  1 SW of Fuego volcano. The model supports a deeper, smaller magma chamber that probably fed the upper chamber; the assigned density (2.36 g/cm 3 ) assumes that it is still molten and above the Curie point temperature. Two dashed lines extend from the positions in the surface in which L3 intersects the rim of the structure down to the vicinity of the intrusion. The location of L0 along this profile is also shown. body of those dimensions was emplaced underneath the area. Subsequently, magma injection probably stopped and the body started cooling down increasing its density and inducing local gravitational collapse. Topographic, morphologic, and gravimetric features tend to support this suggestion.
In Figure 7 the two rim faults intersected by L3 extend from the surface down showing a good correlation with the dimensions of the modeled, cooled magma body. The density assigned to this body is 2.78 g/cm 3 , slightly larger than that associated with the limestone formation enclosing it. The model supports the existence of a deeper, smaller magma chamber of density 2.36 g/cm 3 , similar to the lower chamber of Nevado volcano that probably fed the overlying magma body and that, in the model, it is also assumed to be above the Curie point temperature.
In order to test the hypothesis that erosion is responsible for the existence of La Escondida depression [19] [21] we filled out the collapsed region with volcanoclastic materials of density 2.55 g/cm 3 and obtained a maximum anomaly amplitude increment of only 3 mGals at the center of the structure with the amplitude difference decreasing toward its borders. This implies that such mass deficiency (i.e., the mass removed by erosion) alone cannot account for the observed gravimetric anomaly. As observed in Figure 7 the amplitude of the negative, gravimetric anomaly associated with this structure is of ~10 mGals. To further test the model we assigned a density of 2.75 g/cm 3 (limestone) to the two uppermost layers (volcanoclastics), both of densities 2.55 g/cm 3 ; model calculations showed that the negative anomaly practically disappears under such conditions. This suggests that, most likely, limestone was replaced by volcanoclastics in those layers as it was displaced downwards.

L4 Southern Magma Body
The model of the gravimetric low along L4, SW of La Escondida structure is shown in Figure 8. Given its The density assigned to this body is 2.36 g/cm 3 , the same as that of the lower magma chamber along L3. In this cross-section the intersection with L0 is located in the middle of the magma body; Cretaceous limestone completely encloses it.
To the E a body of density 2.5 g/cm 3 is necessary to adjust the associated gravity low. proximity to the Colima Volcanic Complex it is interpreted and modeled as a magma body. A red, closed line around the intersection of L0 and L4 in Figure 10 shows the projection, on the surface, of the perimeter of this body. We also observe on the surface a 350 m topographic bulge that we interpret as inflation induced by magma injection under this region. Figure 8 shows the model along L4, transverse to L0, for the inferred Southern Magma Body (SMB); in the scale of the model it is difficult to appreciate the bulge, but at larger scales its presence becomes quite clear. The assigned density (2.36 g/cm 3 ) is the same as that of the lower chamber under La Escondida structure; other densities are shown in the model and in Table 1. The magma body is considered to be above the Curie temperature and thus it is modeled with no magnetization. In the model the intrusion is completely enclosed by Cretaceous limestone. To the east of this magma chamber the model requires a smaller body of density similar to the overlying volcanic formations; however, the nature of this body cannot be ascertained with the present information. The magnetic susceptibilities of the associated geological formations are shown in Table 1. As in the previous models, the two intersecting profiles were modeled independently, thus reducing the uncertainty in the geometry and the parameters characterizing it. This magmatic body thus appears as an irregular disk of radius 10 km and 1.5 km average thickness buried at a depth of 1.5 km below sea level. The depth of burial and the thickness are similar to those obtained for La Escondida. In fact La Escondida and the SMB appear along L0 as if they were a continuous body; in the next section we discuss the model along L0 that allows for a better appreciation of this statement. The result obtained by [9] on his magnetic model line extending along the two summits, could have included the rim of the SMB but his depth estimate is about twice the one obtained here.

L0 Base Line
The models along lines L1 through L4 will now be linked by the transverse line L0 that incorporates the densities, magnetic susceptibilities and geometries of those lines in the corresponding intersections. In Figure 9 we show the model associated with line L0. The model is 58 km long and reaches a depth of 15 km crossing lines L1 to L4 in a nearly perpendicular direction. The areas shown in red and magenta are the low-density bodies that we associate with magmatic chambers of either 2.30 or 2.36 g/cm 3 densities; these values are similar to those usually assigned to molten silicic magma [7]. Some magmatic conduits presumably feeding the chambers are also shown; however, owing to their small thicknesses and emplacement depths they are not well constrained. The active conduit to the summit of the Fuego volcano does not reach to the surface since this line (L0) does not go exactly through that point. Magmatic chambers appear at two depth levels; the uppermost level is emplaced in the Cretaceous limestone present in the area at an average depth of 1.3 km below sea level. Those in the lower level appear to be under the bottom of the limestone formation. The densities of the upper volcanic formations vary from 2.42 to 2.45 g/cm 3 . Densities are identified in each formation in this figure and are related to the various volcanic products found in the area. The whole region appears to be underlain by basement rock of densities 2.80 -3.10 g/cm 3 , through which the magma has penetrated to the upper layers. Regarding the fitting of the aeromagnetic anomaly we observe that, as expected, the uppermost layers dominate the magnetic response. The magnetic parameters associated with them were eliminated from the figure in order to avoid cluttering; they are presented in Table 1.
SW of the magma chambers of Nevado and Fuego volcanoes there is another group of chambers at depths of emplacement similar to those of the above volcanoes. They do not have volcanic edifices associated; however, the surface manifestations of their underground presence are a collapse structure discussed previously with respect to La Escondida structure in L3 and a 350 m bulge in the area of the southwestern most chamber. The model corresponds to a sill-type structure that we interpret as another magma body that is in the process of solidification at its northern end (La Escondida) while the southern part (SMB) is interpreted as molten material. Similar intrusive bodies are not uncommon in volcanic fields yielding either positive or negative gravimetric anomalies [7] [24] [25]. In Figure 9 one can also see that the body associated with this anomaly is modeled as an intrusion of density 2.78 g/cm 3 whose upper surface is the horizontal continuation of the neighboring magma body of density 2.35 g/cm 3 . This apparent continuity suggests that these two bodies may have a common feeder, which first injected La Escondida portion and then migrated to the south. As pointed out earlier, La Escondida corresponding magma body probably cooled off increasing its density and inducing collapse, creating the rim observed at the surface.
The model of the CVC along this line suggests a complex, interconnected magmatic system, spreading in the N-S direction. The top of the three magmatic chambers is located ~1.2 km below sea level and their bottoms are located within 3 to 6 km below sea level; these locations and dimensions are fairly similar to those obtained by [9] for the magma chambers of the Fuego and Nevado volcanoes. The magma chamber under Fuego volcano is the only one presenting an active conduit to the surface. Under Nevado volcano the model suggests the presence of at least two, nested, magmatic chambers whose conduits to the surface are now blocked. The assigned densities (2.35 g/cm 3 ) are slightly higher than the one assigned to the Fuego volcano and their magnetizations are considered nil, implying that although they are probably cooling off, their temperatures are still above the Curie point. Figure 10 shows the DEM model of the CVC, lines L0 through L4, and the projections of the modeled bodies on the surface as determined from the 2-D models presented above. Three of them correspond to magma bodies and one to the intrusion under La Escondida structure. Also shown are the fault scarps of the graben-like structure inferred by [19] to reflect the geometry and style of deformation in the edifice interior. The correlation between the sizes of the modeled bodies and the corresponding topographic collapses in the surface confirm the Figure 10. Surface extent of the modeled bodies, molten (red: Fuego and SMB) or solidified (magenta: Nevado, yellow: La Escondida), projected onto the DEM of the region [13]. The modeled lines L0 through L4 and the gravimetric stations and contours are shown for reference and are the same as in Figure 2. Three triangles in the center represent the locations of Los Hijos, Fuego, and Nevado structures (from S to N) respectively. The two additional triangles in the N portion correspond to Cántaro S and Cántaro N. Dark blue traces correspond to the graben structure of Nevado, and light blue traces correspond to collapsed regions, as reported by [19]. The collapsed region of La Escondida is shown in green. A good correlation is observed between the dimensions of the modeled bodies and those of the corresponding surface structures. The magenta polygon defines the surface area under which 3-D inversions are performed.

Surface Projections of Modeled Bodies
close relationship between them as described by the latter authors.
In the case of Nevado there is an excellent correlation between the dimensions of the chamber and the deformations of the structure, including the E-W graben.
Correlation of Fuego volcano and the caldera collapse is also quite close.

Gravimetric and Magnetic Inversions
The 2-D models presented above will now be complemented by two 3-D inversions of the magnetic and gravimetric fields shown in Figure 2 and Figure 3. We use the inversion method described by [26] based in turn in the theoretical considerations of [27]. A 3-D mesh is built under a selected surface area in which a grid of magnetic or gravimetric values is defined. In the present case the selected surface area encompasses the volcanic structures and is delimited by the magenta polygon in Figure 10. Each volume element (voxel) of the mesh is assigned a density or a magnetic susceptibility, depending on the type of inversion performed. Next, the contribution to the total field of each voxel is calculated at an observation point in the surface. The total contribution of the set of voxels is compared to the corresponding measured value on the surface. The process is repeated for each observation point in the surface until the difference between the calculated and measured values is less than a predetermined value. In our case this value was chosen to be 5 percent of the standard deviation. Thus the grid generated by the inversion process does not differ from the measured grid more than 5 percent of the standard deviation.

Gravimetric Inversion
Inversion of the gravimetric data will first be discussed. The box in Figure 11(a) shows the volume in which calculations were made; prisms dimensions were X = 1000, Y = 1000, and Z = 500 m and range in depth from 4200 to −4500 m, thus yielding a blocky appearance to the topography. The density scale is relative to the reduction density (2.67 g/cm 3 ); the low value in the scale (−0.22) corresponds to 2.45 g/cm 3 , while on the opposite end (0.05) it corresponds to 2.72 g/cm 3 . Some of the prisms have been suppressed in order to show a vertical cross-section along L0, complementing the 2-D model along this line (Figure 9). The assigned densities in the 2-D models and the densities obtained in the inversion process are quite similar. In spite of the coarser resolution of the 3-D inversion there are important remarks about the density distribution in the region. One can readily identify the two vertical regions (CVC and SMB) in which the 2-D models located the presence of magma chambers; here they correspond to low-density regions, associated with the magma bodies. Inversion of the gravimetric data will first be discussed. The box in Figure 11(a) shows the volume in which calculations were made; prisms dimensions were X = 1000, Y = 1000, and Z = 500 m and range in depth from 4200 to −4500 m, thus yielding a blocky appearance to the topography. The density scale is relative to the reduction density (2.67 g/cm 3 ); the low value in the scale (−0.22) corresponds to 2.45 g/cm 3 , while on the opposite end (0.05) it corresponds to 2.72 g/cm 3 . Some of the prisms have been suppressed in order to show a vertical cross-section along L0, complementing the 2-D model along this line (Figure 9). The assigned densities in the 2-D models and the densities obtained in the inversion process are quite similar. In spite of the coarser resolution of the 3-D inversion there are important remarks about the density distribution in the region. One can readily identify the two vertical regions (CVC and SMB) in which the 2-D models located the presence of magma chambers; here they correspond to low-density regions, associated with the magma bodies.
Both regions widen with depth; this suggests they may be connected at greater depths. The region corresponding to the CVC locates the lowest density (−0.27 or 2.40 g/cm 3 ) at its bottom (−4500 m) also showing a good correlation with the results of the 2-D models. A similar observation is made at the SMB although the density is slightly higher. Both regions appear to increase density in concentric layers. The CVC and SMB are separated by a region of higher densities (blue tones) that does not quite reach to the surface. Both anomalies start abruptly in the N and then decrease gradually towards the S. The 2-D models showed the upper magma chambers located at burial depths of ~1200 -4300 mbsl. In Figure 11(b) we show a slice of the inverted density distribution between depths of 1100 and 2100 mbsl in the region of the CVC gravity anomaly. The topography derived from the DEM is superposed as a reference. The slice shows lowest densities under Nevado and Fuego edifices and the SMB, as well as higher densities in La Escondida region (E). A SW-trending channel connecting the two low-density zones can also be appreciated as well as some low-density ramifications to the NW in the zone corresponding to the CVC and to the E and W in the SMB. This channel corresponds to that previously discussed in connection with Figure 2. The high correlation between the 2-D models and the gravity inversion results is evident.
A 3-D representation of the lower-density voxel distribution is shown in Figure 11(c). All voxels of densities >2.45 g/cm 3 were clipped off from the 3-D inversion; the remaining voxels are located between 0 and −4500 mbsl. Two low-density regions remain corresponding to the CVC and the SMB after the filtering process. Thus  Figure 2 in the area defined by the magenta polygon shown in Figure 10. The cross-section shown is along L0. Red represents low densities and blue high densities. The scale represents densities relative to the reduction density (2.67 g/cm 3 ). The box delimits the volume in which the voxels used for the calculation are defined. Voxel dimensions are X = 1000, Y = 1000, and Z = 500 m. The Z-range is from 4200 to −4500 m. Low-density regions are associated with the CVC and the SMB; the density corresponding to the latter is somewhat larger than the former. These density distributions spread with depth, and in the upper portions appear to vary concentrically. (b) Density distribution slice between depths of 1100 and 2100 mbsl with DEM superposed (vertical exaggeration 1.5, [13]). The 2-D models located the upper magma chambers of Nevado, Fuego and La Escondida (E) structures between those depths. Hence, density distribution obtained from the inversion process correlates well with densities in the 2-D models. CVC Colima Volcanic Complex, E La Escondida, CN Cántaro N, CS Cántaro S, N north, S south. (c) All voxels with values >2.45 g/cm 3 were clipped off the distribution shown in Figure 11(a); the voxel color convention is the same. Those voxels remaining correspond well with the locations of the magmatic chambers of the CVC and the SMB. The DEM [13] is superposed for correlation purposes. Resolution is too coarse to resolve the independent Nevado and Magma chambers.
we are looking at the two regions likely to contain the magma bodies; both show a dome-like structure and no connection between them is observed at these depths and these densities. The DEM representation, with a vertical exaggeration of 1.5, is superposed for correlation and reference purposes; it is clipped off along a line parallel to L0 and slightly to the W of it, in order to allow for a full view of the low-density bodies underneath. We note that with the voxel size used in this inversion the magma chambers of the Nevado and Fuego volcanoes cannot be resolved; they appear joined in a single, low-density body.

Magnetic Inversion
Based on the aeromagnetic map of Figure 3 an inversion was made similar to the gravimetric inversion presented above, with the same dimensions of the inverted area and the same dimensions of the prisms used. Therefore, the gravimetric and magnetic inversions presented here are directly comparable between them as far as location and resolution is concerned. The physical property yielded by the magnetic inversion is the magnetic susceptibility. Figure 12(a) shows a cross-section of the inversion along L0 similar to that presented in Figure  11(a) for gravimetry. We observe two magnetic susceptibility regions of ~0.002 SI which are interpreted as arising from the geologic materials surrounding the magmatic chambers; a kind of magnetic susceptibility halo. In between these two regions there is a susceptibility low region that seems to correspond to the intrusive associated with La Escondida structure. The susceptibility layers appear to be in concentric layers, similarly to the concentric distributions observed in the corresponding density layers. We infer thus, that the susceptibility distribution also outlines the regions lodging the magmatic chambers. Figure 12(b) exhibits a 3-D representation of the magnetic susceptibility surfaces with values of 0.000 and 0.007 SI extracted from the inverted susceptibility distribution. The former value is pertinent since it corresponds to the susceptibility value assigned to the magma chambers in the 2-D models ( Table 1) yielding information of their location from the magnetic viewpoint. Two cylindrical regions are delimited by the orange surfaces; one corresponds to the region below the edifices of the CVC and the other below the location of the SMB. A third orange region to the NW is marginal in position and will not be discussed here. Another susceptibility surface (magenta) is located within the large cylindrical region of the CVC that shows two hollow regions in the upper part. The positions of these hollow regions correspond quite well with the corresponding surficial positions of Nevado to the N, and Fuego volcano to the S suggesting that this surface may limit the presence of the magma chambers of these edifices. We noted previously that, owing to the low resolution of the gravimetric inversion, one cannot discern between the magma chambers of Nevado and Fuego volcanoes. However, it appears that combining the gravimetric and magnetometric inversions a better definition of the location of the chambers can be made. In order to show the excellent correspondence between both inversions we added to Figure 12(b) the distribution of voxels shown in Figure 11(c), corresponding to the magma chambers of the CVC and the SMB. These voxels appear as a wire mesh representation in order not to obliterate the susceptibility surfaces. As can be seen, the cylindrical surfaces enclose quite well the density voxels interpreted to pertain to the magmatic chambers. These inversions involve different potential fields independent of each other, and correspond to different physical rock properties. For these reasons we submit they fully support our 2-D and 3-D interpretations.

Discussion
As previously noted, various authors have pointed out the anomalously close position of the CVC with respect to the Middle America Trench (MAT) [1] [28] [29]. Although the reason for this proximity still waits for a full explanation, [4] points in a direction that may help clarify it, and also elucidate the problem of the southward-younging activity in the CVC. The latter authors, analyzing seismic horizontal, tomographic slices in the region found that while the Rivera and Cocos plates arrive jointly at the MAT they begin to separate at depths of 140-200 km; furthermore, the surface location of the place at which the Rivera and Cocos plates begin to separate at depth according to their results is slightly to the west of the Central and Northern Colima rifts, which they assume to be the boundary between the plates at this location. The larger dip of the Rivera plate under the Jalisco block, compared to that of the Cocos plate, has also been noted (e.g., [4] [30]), modifying the rather uniform subduction regime existing along this margin prior to the separation of Baja California from mainland Mexico [31]. The Northern Colima rift started to form in the Late Miocene [23]; however, the oldest volcanic activity of the CVC is registered at only 1.7 Ma in El Cántaro volcano [2], indicating that deep tectonic processes responsible for the rift formation took about 4 M years to manifest themselves at the surface in the form of volcanism.  Figure 11(a). The region of the CVC widens with depth while that corresponding to the SMB does not. (b) 3-D representation of two magnetic susceptibility surfaces with values of 0.000 (orange) and 0.0018 SI (magenta) extracted from the inverted susceptibility distribution of Figure 12(a). The former is the susceptibility assigned to the magma bodies in the 2-D models. Note the hollow regions defined by the magenta surface that correspond quite well with the surficial positions of the structures of Nevado to the N and Fuego to the S. The two groups of low-density voxels shown in Figure 11(c) are incorporated in this figure to highlight the excellent correlation between the positions of the magma voxels derived from the density inversion and the limiting magnetic susceptibility surfaces derived from the magnetic inversion. The voxels' wire mesh representation is used in order to avoid surface obliteration.
Reference [5] confirms that the Rivera and Western Cocos plates have rolled back towards the trench and that the Colima rift is intimately related to the tearing between those plates. We can only speculate that at some point in time the gap at depth between the Cocos and Rivera plates became large enough as to allow toroidal flows [4] to take place, some of its products migrating to the surface producing volcanism. If the Rivera plate is rolling back to more vertical positions [4] [32] [33] then the gap that allows toroidal flows between the Rivera and Cocos plates is also growing and displacing southwardly, providing a plausible explanation to the southwardyounging volcanic activity of the CVC. If this is the case, however, there remains to be explained why volcanism ceases in well-developed volcanic structures (e.g., Cántaro and Nevado) when the southward migration occurs, instead of maintaining the magma flows to their respective magma chambers from an ever growing gap between the subducting Rivera and Cocos plates.
At present the volcanic activity of the CVC is located in the Colima volcano while Nevado de Colima is considered an extinct volcano. However, there is evidence of successive stages of activity in the latter edifice as recent as ~250,000 years BP, with the dating (K/Ar) of the Atenquique conglomerate derived in turn from La Calle andesite [14] [15]. This andesite outcrops at the top of Nevado in an area of ~20 km 2 . Reference [15] reports an age of 9370 ± 400 yr B.P. for charcoal debris directly overlying debris avalanche deposits SE of the caldera of the ancestral Colima volcano implying that construction of the ancestral Colima volcano edifice took place before that date. In fact, these authors established a period in which the two modern structures of Nevado II-III and Paleo Fuego-Fuego were active, indicating that their contemporaneous growth poses the problem of the existence of a single reservoir feeding them, or the existence of two independent magmatic chambers after the evolution of a common, primitive volcano. The period of simultaneous activity in both structures may correspond to the transition period in which volcanic activity was migrating from north to south. Our models reach depths of ~20 km along L0 and suggest the existence of independent magma chambers for each volcanic edifice, which appear to be interconnected at depth. Our model thus favors the possibility of an even deeper magma chamber feeding all the more superficial ones. This deeper magma chamber is probably connected with the tectonic process that separates at depth the subducting Rivera and Cocos plates [4]; the plates begin to separate between depths of 140 -200 km according to their horizontal, tomographic slices. However, one must bear in mind that the degraded resolution of their method at depths less than 100 km may obscure shallower phenomena. As pointed out earlier, their surface location of the place at which the Rivera and Cocos plates begin to separate at depth is slightly to the west of the Central and Northern Colima rift.

Conclusions
The concurrent modeling of the gravimetric and magnetic anomalies reduces the ambiguities inherent to these potential field methods and results in an improvement over previous models based on the single inversion of either one; in addition, we have extended the length of the base line L0 with respect to the previous models and added smaller, perpendicular lines over the main structures. We located the two magmatic bodies related to Nevado and Fuego Colima volcanoes along the base line and the two sub perpendicular lines L1 and L2, determining their extensions and depths of burial. The size of the Nevado and Fuego modeled magma chambers correlates well with the sizes of the corresponding collapsed calderas; in addition we have extended the surveyed area around the volcanic structures and found additional gravimetric anomalies to the south of the CVC. Inversion of the observed magnetic and gravimetric fields yielded 3-D susceptibility and density distributions which are concordant with the 2-D models. Of particular interest is the body identified as the Southern Magma Chamber since it is nearer to Colima City (population 250,000); at this point it is not possible to predict the evolution of this magma body either into a new volcanic edifice or another collapse structure.
The observation regarding the southward-younging activity in the CVC can now be extended to its southern limits. The sequence and proximity of Nevado volcano, Fuego volcano, La Escondida structure, and the Southern Magma Chamber favor a common origin for them, and its alignment suggests that the southward migration of volcanism is an active, ongoing process. This is in agreement with the conclusion [1] that Colima's andesites are the product of a dynamic subvolcanic magma system, fed from below by pulses of relatively basic magma. However, the question remains of why an active volcanic structure is abandoned after trench ward emplacement of new magmatism.