Introduction

Study of radiometrically dated lava flows permits understanding the behavior of the Earth magnetic field (EMF) during the geological periods. The paleomagnetic records provide valuable information about the directions and intensity stored in each independent cooling unit. Sediments may provide quite continuous records of magnetic field variation, while lavas, due to the sporadic character of volcanic eruptions, yield rather discontinuous records of geomagnetic field fluctuations. On the other hand, the results obtained from lavas are generally more reliable because of nature and physical principles of the thermoremanent magnetization (TRM) acquisition (^{Prévot et al., 1985}). Many high resolution volcanic records have documented polarity transitions and intervals of constant polarity termed chrons and relatively short duration of 103−104 year’s events or excursions inside. The transitional episodes are generally defined in terms of a deviations in intervals less than 103 years of the Virtual Geomagnetic Pole (VGP) position from the Geomagnetic Axial Dipole (GAD).

The whole Trans Mexican Volcanic Belt (TMVB) represents and excellent target for high standard paleomagnetic studies since it offers more than 3000 Plio-Quaternary lava flows (mainly monogenetic volcanoes) being many of them radiometrically dated using either K-Ar or Ar-Ar systematics. The present investigation is aimed to contribute to the Time Averaged Field (TAF) study for the past My and improve knowledge on the paleosecular variation (PSV) at low latitudes (20°) extending the previous work of ^{Maciel et al. (2009)} from Tancitaro volcano and the surrounding MichoacanGuanajuato Volcanic Field (VFMG).

^{Johnson et al. (2008)} analysis of a new generation of paleomagnetic data provides new insights about the latitudinal dependence of VGP (Virtual Geomagnetic Pole) scatter showing that this relationship is much less comparing to previous studies (^{Tauxe and Kent, 2004}). It appears that the latitude dependence of VGP angular dispersion depends critically on a dataset from the moderate latitudes about 20° north including numerous data from Hawaii and Central Mexico. In other hand, modeling of the aeromagnetic anomalies in the area permits to investigate the subsurface structure, and stratigraphy, in particular about distribution of magnetic sources related to magmatic bodies and volcanic structures.

Geological setting and sampling details

Tancitaro volcano (TV) is an andesitic-dacitic stratovolcano located in the MichoacánGuanajuato volcanic Field (MGVF). The MGVF has an extension of 40,000 km² and is part of the central portion of the TMVB. The MGVF has geographic boundaries 18º45’N and 20º15’N in latitude and 100°25’W and 102°45’W longitude; contains over 2000 small-sized monogenetic volcanoes including basaltic monogenetic cinder cones (^{Hasenaka and Carmichael, 1985}), maars, tuff rings, lava domes and lava flows with hidden vents. Volcanic products are predominant calc-alkaline, but some alkaline and transitional rocks are present. Silica content varies from 47% to 70% for olivine basalt and basalt-andesite rocks (^{Hasenaka, 1994}; ^{Hasenaka et al., 1994}).

The TV is a highest stratovolcano in the MGVF with a height of 3840 mts (^{Ownby et al., 2007}), and has a large amount of dated flows. The activity start ≥793±22 kyr, and its last eruption was around 237±4 kyr, most of these events were dated by the radiometric method Ar⁴⁰/Ar³⁹ (^{Ownby et al., 2007}; ²⁰¹⁰). The present study sampled 8 lava flows to complement previous work made by ^{Maciel et al. (2009)}. These lava flows (Figure 1) have been reported by ^{Ownby et al. (2010)} who present 39 new Ar-Ar dates that complement 26 flows with radiometric age also reported by ^{Ownby et al., (2007)}.

Figure 1 Location of the Tancitaro Volcano in the western sector of the TMVB showing the setting of sites reported in this study (Tan_01-08) and the previously sampled flows by ^{Maciel et al., (2014)} (Tan1-11). 

For the present study, the volcanic flows correspond to those with radiometric age, easy access, and sampled roadside of almost continuous structures of lava flows not being rotated or unaltered outcrops and not sings being near a geological fault which was corroborated with the geological map of the National Geological Service. The eight samples were obtained with a gasoline powered portable drill, oriented with a magnetic compass to obtain several samples along each flow N≥8.

The eight flows sampled cover an age range of ~70 to 957 kyr, and complements previous work by ^{Maciel et al., 2009}, that included 11 flows with an age range from 82 to 612 kyr, since the TV is an important source of unrecorded flows that could provide important paleomagnetic information with aged of a good accuracy.

Magnetic Measurements and Data Analysis

Remanence measurements

All specimens were demagnetized by peak alternating fields which proved to be highly efficient to isolate the characteristic, primary remanence. A Molspin AF-demagnetizer with available alternating fields from 5 to 95 mT was used while magnetization was measured with a JR6A (AGICO) spinner magnetometer with nominal sensitivity ~10−9Am2. The determination of the main magnetization components for each specimen was achieved with the method of principal component analysis (^{Kirschvink, 1980}) and the directions were averaged by unit based on Fisher statistics (^{Fisher, 1953}).

In most of the cases a stable, uni-vectorial characteristic remanent magnetization (ChRM) was observed (Figure 2 a) for Tan1 and Figure 2 d) for Tan7) occasionally accompanied by a negligible viscous overprint (Figure 2 b) for Tan 2 and Tan5) easily removed after the first step of demagnetization (Figure 2). Most of the specimens were completely demagnetized until 100 mT. The site mean directions are quite precisely determined since all a95 are found less than 10°. For one out to eight studied lava flows (Tan4) not paleodirection are determined due the erratic and unstable behavior during the magnetic cleaning of specimens. The other seven volcanic flows gave stable paleodirections, and were divided in two groups, the first group are all the directions that belong to the chron of Brunhes (5 of the 7 paleodirections) and the second group are the paleodirections that belong to the chrons of Matuyama. The directions that belong to the chron of Matuyama had the peculiarity to correspond to the same event of transition due the normal polarity and the radiometric age assigned to these flows by ^{Ownby et al. (2010)}. The transitional event known as Jaramillo, was first recognized by ^{Doell and Dalrymple (1966)}.

Figure 2 Representative examples of demagnetization experiments using peak alternative fields up to 95 mT, with the orthogonal vector plots, (Zijderveld diagrams) a) For the sample 94T004A and the site Tan1 b) Sample 94T015A and site Tan_02, c) Sample 94T047A corresponding to Tan_05 and d)Sample 94T060A for the site Tan_07.  

Rock magnetism

The acquisition of thermomagnetic curves for representative samples are reported on Figure 3. Curie temperatures were estimated using the differential method of ^{Tauxe (1998)} for the analysis of the Magnetic vs Temperature curves. These analyses evidenced the low temperature phase in the heating process (red line Figure 3) Tc 492°C for sample 94T002A corresponding to site Tan_01 (Figure 3a) with two phases in the heating process (red line) that may indicates the presence of Titanomagnetite with medium Ti content. Using the method of ^{Moskowitz (1981)}, it results that the highest Curie temperatures are 591°C for 94T016A corresponding to site Tan_02 (Figure 3b), and 511°C for 94T028A corresponding to Tan_03 (Figure 3c). In some cases the heating and cooling curves are not perfectly reversible (Tan_01). This may be due to the high to moderate mineralogical alteration at high temperatures occurred during the laboratory heating.

Figure 3 A summary of rock-magnetic experiments for the most representative samples: Temperature dependence of magnetization (magnetization ), the red curve represent the heating process and the blue the cooling curve. Hysteresis loops for an induced magnetic field in blue (magnetization). Isothermal remanence magnetization acquisition curves (IRM, )obtained with Variable Field Translation Balance with the respective Day Plot (^{Dunlop, 2002}), to estimate the domain state of magnetic carriers and the relation of the ratios of the hysteresis parameters for the remained samples. 

The hysteresis loops were analyzed using a RockMagAnalyzer1.0 software (^{Leonhard, 2006}). Near to the origin, there is no evidence of wasp-waisted behavior (^{Tauxe, 1998}), which reflect restricted coercitivity ranges. When judging the ratios obtained from the hysteresis curves, samples fall in the pseudosingle domain PSD field (Figure 3). Isothermal remanence acquisition curves are sensitive to the magnetic mineralogy, concentration and grain size properties. Almost all samples are saturated at about 300 mT applied magnetic field, which indicate the presence of ferromagnetic phase with moderate coercivity as may be expected from magnetite and titanomagnetite grains (^{Tauxe, 1998}). The lack of total symmetry on the hysteresis plots is generally attributed to the instrumental errors except for the study performed by ^{Chandra et al. (2012)} where the exchange bias seem to be the principal cause.

Main Results and Discussion

All lava flows associated to the TV yielded a normal polarity magnetization. The paleodirections were divided in two groups, five of them (Tan_01, Tan_02, and Tan_06-08) correspond to Brunhes Chron, and the other two (Tan_03 and Tan_05) belong to Matuyama Chron.

The mean paleodirection obtained from the seven flows is, with an uncertainty α₉₅=8.21º and K=13.89, which correspond to the pole position of Plat=83.06º, Plong =294.86º, A95=14.49º. These paleodirections are combined with those obtained by ^{Maciel et al. (2009)} that reported 11 normal polarity directions within normal Brunhes Chron (Table 1, Figures 4 and 5). For the global mean calculation purpose we selected sites with N > 4 samples per site and α₉₅ ≤10º. Moreover, transitional polarity data are rejected as common in studies of paleosecular variation (PSV) and time average field (TAF) for an age less than 5 Myr (^{Johnson et al., 2008}, Ruiz-Martínez et al., 2010) to obtain a mean calculation.

The mean paleomagnetic direction obtained in this study is I_m=37.46º, D_m=358.39º, α_95m =6.4º, K =34.15 N=15 and (Figure 5a). The corresponding paleomagnetic pole position is P_lat=86.57º, P_longm=246.34º, A_95m=6.21º and K_m=33.87 (Figure 5b).

Table 1 Summary of all the data used in this article. Site: are the label of the sampled site; Dec and Inc: are the paleodirection treated in this study; k and a95: parameters of quality within the cone of 95%; N: Is the number of samples used to obtain the paleodirection; latitude and longitude are the coordinates for each one of the sampled sites. Age: represent the age obtained by the Ar-Ar method with its respective dispersion; VGP: are the virtual geomagnetic pole for latitude and longitude; Ref: Differentiate the paleo-direction of this study (T.S) and the paleodirection obtained by ^{Maciel et al. (2009)}. 

Figure 4 The corresponding VGP of the paleodirection with the mean VGP, and compared with the poles of reference for the last 5 Myr. 

Figure 5 a) Projection of the paleodirection for the Tancitaro volcano with the respective mean paleodirection compared with the direction of reference for the Craton of North America provided by ^{Besse and Courtillot (2002)} and ^{Torsvik et al., 2012}. b) The corresponding VGP of the paleodirection with the mean VGP, and compared with the poles of reference for the las 5 Myr. 

The mean paleodirection obtained was compared with the expected direction derived from the reference pole positions for North American Craton proposed by ^{Besse and Courtillot (2002)} and ^{Torsvik et al., (2012)} yielding D_B =356.1, I_B =33.1 and D_T =358.5, I_T =34.2 respectively. The use of North American reference poles is a common practise in almost all previous paleomagnetic surveys (^{García-Ruiz, et al., 2016}, ²⁰¹⁷; ^{Maciel et al., 2009}, ²⁰¹⁴; Ruiz-Martínez et al., 2010.) in order to avoid a kind of bias due to the neotectonic acivity. For instance, the San Juanico Buenavista Fault located just beneath Tancitaro Volcano may potentially produce some minor tectonic rotations.

The mean direction agree with the expected directions within the uncertainty and this may be corroborated with the angle of deviation δ_B =4.54 (with respect to ^{Besse and Courtillot (2002)} reference pole) and δ_T =3.34 (with respect of ^{Torsvik et al. (2012)} poles). For each paleomagnetic study, it may result useful to calculate the Flattening and Rotation (^{Butler, 1991}) parameters together with their corresponding statistics. These analysis indicate the absence of any important tectonic movement since the values are relatively low R_B =1.66±7.22, F_B =−4.44±6.33 (for ^{Besse and Courtillot, 2002}) and R_T =−0.73±7.22, F_T =−3.34±6.33 (for ^{Torsvik, 2012}). These values also agree to the study carried out by Ruiz-Martínez et al. (2010) at different areas within the TMVB, specially for the western and central sectors, where similar small values are obtained ( ⎪R⎪<5º).

Most paleomagnetic studies make the implicit assumption that when averaged over some time interval, the paleomagnetic directions are close to the geocentric axial dipole (GAD) for at least the last 5 Myr, as globally averaged paleopoles should reflect mainly the GAD field (^{Tauxe and Kent, 2004}, ^{Carlut and Courtillot, 1998}). Under this assumption, we developed an analysis of TAF and a statistical comparison with the GAD. The mean direction of GAD was calculated for a mean latitude λ=19º16.08’ obtaining D_GAD = 0º, I_GAD = 34.96º, yielding a difference of Δl=2.48º and ΔD=2.32º. These values are even lower than those found by ^{Maciel et al. (2009)} in previous studies and reasonably close to the GAD directions. The inclination anomaly found in this study is also compared with the value reported by ^{Lawrence and Constable (2006)}, for similar latitudes λ=±20º as Hawaii, Mexico, South Pacific and Reunion. All these regions yielded statistically indistinguishable values (Figure 6).

PSV analyses are very important to estimate the paleofield location for any time known and when a multiple site sampling is collected the influence of the error is reduced. The PSV analysis was developed complementing the directions presented in this study with the directions of ^{Maciel et al. (2009)} and the use of , to obtain S_W=14.09, the standard deviation as S²_p = which result Sp =14.33, the root mean square of the angular deviation of VGP about the geographic axis result S_B =14.33 with respective lower confidence limit of 13.6 and upper confidence limit 15.61, with the use of the ^{Cox (1970)} method.

Figure 6 Curve of the best fit with respect of the anomalies of the inclination for several latitudes, with the average of the inclinations anomalies (blue dots) and the anomaly inclination of Tancitaro (Black dot). 

These results were compared with the model G proposed by ^{McElhinny and McFadden (1997} generated by global paleomagnetic data set of lava flows exclusively. This database was reanalyzed by ^{Johnson et al. (2008)} using the bootstrap method under the main assumption that the paleomagnetic field closely approximates with the GAD for the last 5 Myr (see also a simplified statistical model Tk03.GAD by ^{Tauxe and Kent (2004)}. Our results show a good agreement with the model G proposed by ^{McElhinny and McFadden (1997)} (Figure 7) for a latitude of λ=20º and a dispersion of S_M=14 and it is close to the uncertainty limits of the dispersion for the of ^{Johnson et al. (2008)} model G, with a value of S_J =14.8, but when compared with the Tk03. GAD model some disagreement is observed since this model predicts a low value of S_T =12.9 with respect to the λ=20º latitude.

For the studies developed around of the TMVB, a dispersion S=15.06±0.6º and S =14.9^14.8 was obtained for studies made by ^{Lawrence and Constable (2006)} and Ruiz-Martínez et al. (2010) respectively, with dispersion a little bit higher, due the biased to Brunhes Chron.

Figure 7 Scatter of the VGP as a function of the latitude for the last 5 Myr, compared with the curves of reference provided by de Model G proposed by ^{McElhinny and McFadden (1997)} and ^{Johnson et al. (2008)}, also with the Tk03.GAD model of ^{Tauxe and Kent (2004)}. 

Two paleodirections obtained on sites Tan-3 and Tan-5 yielded normal polarity paleodirections within reverse polarity Matuyama Chron. They may be erupted during the worldwide observable Jaramillo event (^{Doell and Dalrymple, 1966}; ^{Mankinen and Dalryple, 1979}; ^{Laj and Channel, 2007}).

Aeromagnetic data of Tancitaro volcano area was obtained by the Mexican Geological Survey (SGM) in 1999, using the following equipment: Islander Airplane BN2-A21; Geometrics magnetometer G-822A optically pumped cesium vapor with sensitivity of 0.001 nT; data acquisition system, Picodas P-101 AG; Base station magnetometer, GEM Systems GSM-19 Overhauser, sensitivity 0.01 nT; Radar altimeter, Sperry RT-220 Navigation system, GPS Ashtech CG24GPS + Glonnas, 16m. Flight parameters were: contour flights at a height of 300 m; direction of flight lines North-South; distance between flight lines 1,000 m; distance between control lines (East-West) 10,000 m; electronic navigation (GPS). The data was initially processed by a digital compilation for correcting the movement plane (magnetic compensation), and daytime drift. From total magnetic field data the International Geomagnetic Reference Field (^{IGRF, 1995}) was substracted obtaining the Residual Magnetic Field (RMF, Figure 8)

TMF - ^{IGRF (1995)} = RMF. The data were leveled using the control lines and microlevelling. All the processes described above were carried out by the SGM.

From the digital aeromagnetic information of the SGM, we plotted the RMF using the Inclination and Declination data means, in order to calculate the Reduced to the Magnetic Pole (RMP, Figure 9) (^{Baranov and Naudy, 1964}). Based on the RMP, the map of the Derivative in the Z direction was calculated (1 order, Figure 10) (^{Henderson and Ziets, 1949}) and the Upward continuations (^{Henderson, 1970},) among others.

Figure 8 Map showing the calculated Residual Magnetic Field (RMF) of the Tancitaro volcano (TV). The RMF map represents the magnetic field strength with isovalue curves in nT and color range, indicating with red the magnetic highs and with blue the magnetic lows. Note the position of the of the TV with respect to the magnetic anomaly. 

Figure 9 Map of the Reduced to Magnetic Pole (RMP) and Upward Continuation 500 m. The RMP map represents how the study area would be observed if it was located in the north pole, were Inclination is =90º and Declination = 0º. Note the position of the TV with respect to the magnetic anomaly. In this RMP the anomalies move to the north.  

The location of magnetism sources associated with magmatic chambers is based on the analysis of the contour aeromagnetic map (^{López-Loera, 2002}). This analysis consists the identification of normal bipolar magnetic anomaly near or under the volcanic structure. Once the anomaly is identified, the ^{Baranov and Naudy (1964)} algorithm (reduced to magnetic pole) is used in order to locate precisely the anomaly within the studied body. The so called Fast Fourier Transform spectral analysis (^{Blakely, 1995}) was used to estimate the depth of the source. This methodology of identification and association of aeromagnetic anomalies with magmatic chambers is based on the fact that the process of cooling of volcanic bodies start from the outside. In this way a crust is formed around the magmatic chamber that allows its identification.

The information processed was configured obtaining the Residual Magnetic Field (RMF, Figure 8). The mean values for the Magnitude in the aeromagnetic survey data (06/1999) was 41,725 nT, the Inclination 46°31´ and declination 7°45´, which means that the RMF anomaly was displaced from the source associated with them. It is a common practice in the processing of aeromagnetic data to use a mathematical algorithm (^{Baranov and Naudy, 1964}) whose application allows to locate the area with respect to the north pole, where the Inclination is 90° and the Declination is 0° and therefore the magnetism sources will be located below the magnetic anomaly

Contour map of the RMP shows the existence of zones with different magnetization intensities (different colors in Figure 9). This enables to group them into what is called aeromagnetic domains (AMD), López-Loera, 200 TESIS QUITAR2) that can be associated with rocks having similar magnetic susceptibility values. This means that each AMD is associated with a different geological unit. At a regional scale, the entire study area correlates with one AMD, mainly associated with volcanic rocks. At more local scale three aeromagnetic subdomains (AMSD) can be differentiated as described below.

AMSD I. It is located in almost all the study area, it covers the entire central portion and most of the northern side, it is undefined to the E, W and N. It is characterized by the existence of the Tancítaro and Paricutin volcanoes (Figure 9a). It has two large aeromagnetic anomalies (8.8 km x 5 km) elongated in the E-W direction.

Figure 9a Map of the Reduced to Magnetic Pole (RMP) and Upward Continuation 500 m and aeromagnetic subdomains AMSD (I, II and III) of the Tancitaro volcano (TV) and surrounding areas. The AMSD represents a zone with similar magnetic susceptibility. 

The anomalies show intensity values of magnetization from 563 to -271.9 nT. One of these anomalies is correlated with the Tancítaro volcano, whose crater is located 1.4 km SW of the center of the aeromagnetic anomaly. In the contour map of the Residual Magnetic Field, this geological structure shows a normal aeromagnetic anomaly with a 5,450 m of polar distance. According to an analysis of Radial Average Spectrum (Figure 10) the host rock of the magma chamber associated with the Tancítaro volcano is located at depths nearly 2 km., and with the method of average width of the anomaly (Figure 11) it is estimated that the host rock of the magmatic chamber is between 2.6 and 3 km. The Tancitaro volcano shows alignments in all direction, which are associated with weakness zones correlated with faults and/or fractures. Topographically this AMSD I, shows altitudes between 2,031 and 3,792 m, with an average of 2,832 m. Its geological association is with rocks of high magnetic susceptibility such as basalts.

Figure 10 The profile represents the Radially Averaged Power Spectrum of the study area. 

Figure 11 The profile observed with NS direction in the map of the RMP & Upward Continuation 500 m; the same as shown in the lower side of the map and that served as base to estimate the depth of the source of the magnetic anomaly with the average width method. 

AMSD II. It is located across the southern part of the study area and is mainly characterized by a series of aeromagnetic anomalies associated with magnetic minima, showing values of magnetization intensity from -96 to -650 nT (Figure 1, Figure 11). It has an “E” lying form and is undefined to the S, E and W. Topographically it shows altitudes from 1,140 to 2,141m, with an average of 1,785 m. It correlates with rocks having a low magnetic susceptibility, such as volcanic breccia and pyroclasts.

AMSD III. This subdomain is located towards the NW portion of the studied area, and is undefined to N (Figure 12). It has a form tending to a half circle that shows, in its central portion, an anomaly associated with a magnetic high surrounded by magnetic lows. The values of magnetization intensity are between 197 and -257 nT. Topographically it has altitudes of 1,475 to 2,135 m with an average of 1.780 m. Geologically it correlates with rocks with mean magnetic susceptibilities such as andesitic rocks.

Figure 12 First derivative in Z direction and upward continuation 500 m map. This map shows the alignments (white color lines) associated with the anomaly of Tancitaro volcano. This alignments are correlated with faults and/or fractures. 

In summary, the aerial magnetometry on the map of Reduce to Magnetic Pole (RMP) shows the existence of three aeromagnetic subdomains all associated with volcanic rocks. The aeromagnetic anomaly associated with the Tancitaro volcano on the map of RMP is shifted 1.4 km to 50° NE and has an area of 5 km (N-S) x 8.8 km (E-W) and a magnetization intensity of 710 nT. The depth of the magmatic chamber correlated with the Tancitaro volcano is interpreted between 3.2 km and 5 km and its shape is elongated in the E-W direction. It is limited in all directions by alignments correlated with faults and/or fractures.

Conclusions

A detailed rock-magnetic and paleomagnetic study was carried out on the lava flows associated to the Tancitato volcano in order to contribute both to the new generation Time Averaged Paleomagnetic Field database and to estimate the latitudinal dependence of the paleosecular variation trough the virtual geomagnetic pole scatter. The combined dataset incorporating previously reported paleodirections from the same area offer a detailed record of the Earth’s Magnetic Field fluctuation for the last 1 My.

Mean paleodirections was found reasonably close to the Geomagnetic Axial Dipole directions and statistically indistinguishable from the expected directions derived from the stable North America which attest that there is no major tectonic deformations occurred during the last 1 My according to the previous studies. All samples yielded well defined normal polarity magnetization. Two flows are correlated to the Jaramillo polarity event, which provide a useful marker for the volcanic activity in the MGVF. Combined aeromagnetic survey shows the existence of three magnetic subdomains all associated with volcanic rocks. The aeromagnetic anomaly associated with the Tancitaro volcano on the map of RMP is shifted 1.4km to the Northeast 50º and has an area of 5km (North-South) 8.8km (Est -West) and shows a magnetization intensity of 710 nT. The depth of the magmatic chamber correlated with the Tancitaro volcano is interpreted between 3.2km and 5km and the shape is elongated in the East-West direction. It is limited in all directions by alignments correlated with faults and/or fractures, but this could be inside the flows or formed after the events of interest in the geological maps of the Geological National Service, in the paleomagnetic results this relationship was not observed that could be reflected not being near the direction of reference or the GAD.