Zircon (U-Th)/He Closure Temperature Lower Than Apatite Thermochronometric Systems: Reconciliation of a Paradox

: Here, we present seven new zircon (U-Th)/He (ZHe) ages and three new zircon ﬁssion track (ZFT) ages analyzed from an age-elevation proﬁle (Machu Picchu, Peru). ZFT data present ages older than those obtained with other thermochronological data, whereas the ZHe data interestingly present ages similar to those obtained with apatite (U-Th)/He (AHe). It has been proposed that He retention in zircon is linked to the damage dose, with an evolution of the closure temperature from low values associated with a low α -dose (<10 16 α /g), subsequently increasing before decreasing again at a very high α -dose (>10 18 α /g). Studies have focused on He diffusion behavior at high α -dose, but little is known at low doses. We propose that the ZHe closure temperature at α -dose ranging from 6 × 10 15 to 4 × 10 16 α /g is in the range of ~60–80 ◦ C. This value is lower than that proposed in the current damage model ZRDAAM and demonstrates that the ZHe and AHe methods could have similar closure temperatures at low α -dose (i.e., similar ages). These new data strengthen our previous geological conclusions and even highlight a cooling rate approximately twice as important as that deduced from AHe and apatite ﬁssion track data alone at Machu Picchu.


Introduction
To quantitatively unravel the tectonics and/or relief evolution of a given region, lowtemperature thermochronology methods such as (U-Th)/He and fission track dating of apatite (AHe and AFT, respectively) or zircon (ZHe and ZFT, respectively) are often used together e.g., [1,2]. The ZFT and ZHe methods are generally known to record higher temperatures or deeper processes than the AFT and AHe methods [3,4] because of their higher closure temperatures [5][6][7][8]. Additionally, for a given mineral, (U-Th)/He thermochronometry is generally considered more sensitive to lower temperatures than fission track thermochronology [1,5]. Today, these techniques are routinely applied, and numerous studies are published each year for exhumation quantification purposes. However, for over a decade, methodological studies have highlighted that He diffusion in apatites and zircons is strongly dependent on the radiation damage dose, increasing the range of the closure temperature, as important age dispersion with positive AHe-age, and positive and negative ZHe-ages correlations are now often observed e.g., [8,9] for zircons, e.g., [10][11][12][13][14] for apatites.
Regarding the ZHe method, the He diffusion behavior in zircon is still debated. Based on age-effective uranium concentration (eU) correlations, [15,16] proposed that radiation damage produced during U and Th decay influences He retention and loss in zircons, damage produced during U and Th decay influences He retention and loss in zircons, and proposed the ZRDAAM algorithm to model damage production, annealing and He diffusion in damaged zircons similar to RDAAM for apatite [13]. In this model, damage first traps He until a threshold at which point damages coalesce, leading to damage connectivity and the creation of fast He diffusion pathways [15]. In addition, several studies have discussed damage model parameters such as the threshold value e.g., [16,17] or fail to numerically reproduce highly damaged zircons with this classic He kinetic model as the minerals appear to be more He retentive than predicted [10].
On the other hand, [18] proposed a trapping model and a nonlinear relationship between the closure temperature and the increase in α-damage, illustrating the importance of the α-dose for determining the closure temperature. The α-dose corresponds to the total radiation damage that had accumulated in the crystal lattice, mostly due to α recoil. It depends on both eU and the time since when the mineral began to accumulate damage. This model explains not only zircon data presenting age-eU and age-diffusion domain (ESR) correlations, but also data with no such correlation. In the range of geologically possible α-doses, this model predicts closure temperatures consistent with settings within the middle to high range of α-doses (10 17 to 5 × 10 18 α/g), but there do not exist enough available data in low α-dose settings to confirm the model at any existing α-dose. The author of [18] only indirectly suggests a closure temperature for a low α-dose range (<1 × 10 15 α/g) of 60-100 °C with the volcanic zircon data of [19].
In this article, we propose to provide new data with low α-doses to fill the current knowledge gap and improve the above-mentioned model with direct observations. Based on new ZHe, ZFT, and published AHe and AFT data from a Machu Picchu (Peru) vertical profile [20] (Figure 1), we show that in the case of low α-doses, ZHe closure temperatures are closer to AHe closure temperatures than previously proposed, and thus present younger ages than AFT and sometimes even younger than AHe, as predicted by the [18] model. This article does not specifically focus on the geological interpretation of these data but briefly discusses the implications of our results for the geological evolution of the region where samples were collected.

Geological Settings
The Machu Picchu age-elevation profile is located at the center of an Andean morphotectonic peculiarity: The Abancay Deflection [21]. The Abancay Deflection tectonically delimits the Bolivian Orocline to the south (Eocene-early Miocene rotation of up to~65 • ) and the straight and narrow Andes to the north in Peru, Figure 1; [20][21][22]. It partly lies in the Eastern Cordillera (northern part of the Abancay Deflection), where numerous Permo-Triassic batholiths emplaced into Paleozoic metasedimentary rocks of the Marañon complex [23]. Among them, the granitic Machu Picchu Batholith we sampled ( Figure 1) emplaced at 222 ± 7 Ma in the core of the Abancay Deflection [24]. The Eastern Cordillera shows high elevation and relief. The author of [25] proposed that it has been a long-lived structural high because of the absence of a Meso-Cenozoic sedimentary cover. Young AFT and AHe ages were obtained along a quasi-vertical transect of the Machu Picchu Batholith; however, evidence of an unexpected rapid and recent exhumation at a rate of 1.2 km/Myr, initiated at~5 Ma was observed [26,27]. This exhumation phase is also evidenced in another vertical profile 30 km farther east in the Ocobamba valley, still in the core of the Abancay Deflection [27]. In the southern part of the Abancay Deflection, the Altiplano domain is tectonically decoupled from the Eastern Cordillera along the regional crustal-scale Apurimac fault system, Figure 1; [26,27]. In the Altiplano, Eocene plutons   [28,29] emplaced into Meso-Cenozoic sediments [30]. Unlike the core of the Abancay Deflection (Eastern Cordillera), the Altiplano did not experience recent exhumation acceleration, but rather slow and constant exhumation (~0.2 km/Myr) since 40 Ma [27,31]. Because of the differential exhumation pattern, with higher exhumation rates identified at the core of the Abancay Deflection, the active tectonics and curved fault patterns, it has been proposed that the Abancay Deflection should be a tectonic syntaxis comparable to those described in the Himalayan system, for instance [27].

Materials and Methods
We collected seven samples along a 1.9 km quasivertical profile in the Machu Picchu batholith (Table 1). We performed ZHe and ZFT dating from the same samples presented in [20]. In the field, we collected the freshest in situ rocks to avoid sampling near traces of fluid circulation. The samples were crushed and sieved to extract the 100-160 µm fractions in the Géode laboratory (Lyon, France). Zircon crystals were consequently concentrated using standard magnetic and heavy-liquid separation techniques at the GTC laboratory (ISTerre, Grenoble, France).
The zircons were processed at the Dalhousie Noble Gas Extraction Laboratory (Halifax, Canada) for (U-Th)/He dating. They were analyzed following the methods described in [4,5,32]. In parallel, 2 Fish Canyon Tuff standards, 28.48 ± 0.06 Ma; [33] were also analyzed (zFCT-61 and zFCT62 in Table 1). For each sample, 3 single zircon aliquots were run with transparent euhedral grain radii higher than 70 µm, possibly without inclusions and/or fractures. After measurement of their dimensions for α-correction [34] and imaging, each grain was packed into a Nb foil envelope. 4 He was then extracted from each pack in an in-house built He extraction line with successive 15-min-heatings under a focused beam of a 45 W diode laser (1250 • C), until 4 He yields were under 1% of total. After adding a known amount of purified 3 He spike, 3 He/ 4 He ratios were measured with a Pfiffer Vaccuum Prisma quadrupole mass spectrometer. Typical errors are in range of 1.5-2% (1 σ). Samples were analyzed in groups of 36. In each group, 2 Fish Canyon Tuff (FTC) zircon standards were included to ensure accuracy, reproducibility, and reliability of the data. After He extraction, zircons were dissolved in high-pressure dissolution vessels with concentrated HF and HNO 3 at 200 • C for 96 h. Prior to dissolution, samples were spiked with mixed 235 U, 230 Th, and 149 Sm spikes. Isotopic ratios were measured with iCAP Q inductively coupled plasma mass spectrometry (ICP-MS). Additional blank analyses controlled the analytical accuracy. The raw data were reduced using an Helios software package (R. Kislitsyn and S. Stockli). To test different scenarios, we computed the α-dose (Dα) [35] for each zircon dated with ZHe methodology according to the first equation presented in [36]. For ZFT dating at the GTC laboratory, zircon crystals were mounted in PFA Teflon ® and polished [37]. Spontaneous fission tracks were revealed by etching the grain mounts with a NaOH:KOH melt at 228 • C in a covered Teflon dish heated by a laboratory oven for 22-37 h. The criterion for stopping etching was when well-etched tracks were revealed in the majority of the grains in the grain mount. Zircon etching time should be longer [38,39] for low-damaged crystals; however, given the relatively high U-concentrations of 500-800 ppm (Table 2), it is possible to reveal countable tracks in zircons with 10 Ma cooling ages even with intermediate etch times of 15-30 h [37]. Using the external detector method, the samples were irradiated at the FRM II reactor (Garching, Germany). Fission tracks were counted dry on internal grain surfaces at the GTC platform with an Olympus BH2 microscope at 1250x magnification. ZFT central ages were calculated with RadialPlotter software [40], using a ζ-value of 131.49 +/− 5.4 (M. Bernet) for the IRMM-541 uranium dosimeter glass (50 µg/g U).

Results
The six samples collected along the Inca trail (AB-12-64 to AB-17-69) analyzed with ZHe yield mean ages ranging from 2.4 ± 0.3 to 3.2 ± 0.9 Ma, giving a self-consistent age-elevation trend (Table 1; Figure 2). The aliquots of each sample are also consistent themselves. In contrast, the lowest and 6-km-laterally offset sample (AB-17-20) is less consistent: 1 aliquot is young (1.5 ± 0.1 Ma) and fit with the age-elevation trend, whereas the two other aliquots present consistent old ages of~25 Ma ( Figure 2). In addition, the ZHe data do not show any age-effective uranium (eU) content or age-equivalent sphere radius (ESR) correlation ( Figure 3).   [20]; this study. The two upper AFT ages (pink squares) outside of the age-elevation trend were revised for this study.

Figure 2.
Age elevation of the Machu Picchu profile with AHe, ZHe, AFT, and ZFT data [20]; this study. The two upper AFT ages (pink squares) outside of the age-elevation trend were revised for this study.  [20]; thi study. The two upper AFT ages (pink squares) outside of the age-elevation trend were revised fo this study.  (Table 2), are consistent with the other thermochronometers following the same age-ele vation trend (Figure 2; supplementary materials figures S3 to S5).
These new ZHe and ZFT data are much younger and thus apparently incompatibl with the two highest AFT ages (AB-17-68 and AB-17-69 from [20]) in the profile. We thu revised the counting of these two AFT samples (  These new ZHe and ZFT data are much younger and thus apparently incompatible with the two highest AFT ages (AB-17-68 and AB-17-69 from [20]) in the profile. We thus revised the counting of these two AFT samples (Figure 2; Supplementary Materials Figures S1 and S2). In fact, both of the highest AFT samples are difficult to date because of poor apatite quality (fractures, inclusions), as mentioned in [20]; thus, these two samples yield unreliable AFT ages. After revision, the highest sample (AB-17-69) gives an AFT age more compatible with the zircon data (3.1 +3.2 −1.6 Ma) but still not reliable because of the few grains (10) dated and because of the low apatite quality and U zonation. The second highest sample (AB-17-68) with a central age of 18.5 +5.1 −4.0 Ma remains older than the other samples in the profile. Here, also, the very low counts question the validity of this age (Table 3). For these reasons, we will not base the following discussion on the highest fission track data in the profile.

Revisiting the ZHe Temperature Sensitivity
The new ZHe and ZFT age-elevation trends are consistent with the ones obtained with AHe and AFT data by indicating relatively recent rapid cooling, and curiously, the ZHe data are younger or equivalent than usually lower closure temperature thermochronometers such as AHe and AFT (Figure 2). Given the low eU and Dpar values observed in our AHe and AFT data, not abnormally high AHe and AFT closure temperatures should be expected. The ZFT data are consistent with previous data (AHe and AFT; Figure 2) and confirm previous geological interpretations [20,27]. Whereas all ZHe ages, at the exception of the lowest and more distant sample (AB- [17][18][19][20], are younger than AFT and AHe dates from the same samples ( Figure 2). This suggests that those samples present a ZHe closure temperature of~80 • C similar to the AHe system, which is lower than the classically expected, 100-200 • C [1]. One can argue that an analytic issue occurred during the analysis. However, the two Fish Canyon Tuff zircons standards analysis (zFCT-6x in Table 1) give acceptable results (30.7 ± 2.3 Ma and 29.8 ± 2.2 Ma for a 28.48 ± 0.06 Ma age reference [33]), ruling out any strong analytical bias.
Similar observations have been made by [15], who developed the ZRDAAM model, based on correlations between ZHe ages and effective uranium content. They proposed that low closure temperatures are due to middle-low (~10 16 α/g with a Tc =~120 • C) or high annealing damage (>>10 18 α/g). Using the software developed by [43] to implement the ZRDAAM model, we calculated the expected ZHe, AFT, and AHe ages and closure temperatures in function of cooling rates (Figure 4). Modelling yields observed ZHe ages (same or younger than the AFT ages) for cooling rates between 30 and 50 • C/Ma, and indicates that ZHe closure temperature for our samples becomes lower than for the AFT at cooling rates of about 40 • C/Ma. This suggests that the Pliocene cooling might have been faster than previously modelled [20]. at cooling rates of about 40 °C/Ma. This suggests that the Pliocene cooling might hav been faster than previously modelled [20]. To further test the ZRDAAM model [15], we used the HeFTy model [44] in forwar mode in an attempt to reproduce the AHe, AFT, ZHe, and ZFT ages we obtained. We fe the model with a time-temperature history compatible with the one proposed by [20,27 Modelling reveals that ZHe ages are younger than ZFT ages and close to AHe and AFT To further test the ZRDAAM model [15], we used the HeFTy model [44] in forward mode in an attempt to reproduce the AHe, AFT, ZHe, and ZFT ages we obtained. We fed the model with a time-temperature history compatible with the one proposed by [20,27]. Modelling reveals that ZHe ages are younger than ZFT ages and close to AHe and AFT ages, but still older than AHe ages ( Figure 5). The current ZHe model is thus not sufficiently accurate for very low α-doses and cannot robustly explain the observations. Ref. [18] proposed a model similar to the ZRDAAM model in terms of the Tc/α-dose relationship independent of any age/eU or age/ESR correlation or any annealing damage effects. The latter model extends the α-dose range to the low values and indirectly proposes a lower ZHe closure temperature than the ZRDAAM model for very low α-doses, based on volcanic zircon data previously published [19].
In batholiths that cool rapidly, the α-dose upper limit could be approximated by its crystallization age. However, here, the high α-dose could not be geologically explained by the age of the sampled pluton. The granitic Machu Picchu Batholith emplaced at 222 ± 7 Ma [24]. The corresponding α-dose computed for each dated zircon with this emplacement age is between 2.5 × 10 17 and~1 × 10 18 α/g if we assume that all produced α-damage is preserved in the zircon crystals (Table 4). Following the trapping model [18], this value is still too low to allow the ZHe closure temperature to be lower than 80 • C. However, it is interesting to note that the ZFT ages produced on the same samples present young ages <7 Ma, indicating that all fission tracks produced since the crystallization have been annealed. The relatively low observed spontaneous track densities <1.4 × 10 6 tracks/cm 2 despite the relatively high U concentration (600-800 µg/g) of the analyzed zircons indicates that the accumulation of α-radiation damage may not be significant (<4 × 10 16 α/g) in the Machu Picchu batholith zircons (Table 2; Figure 6; [45]), because of relatively high temperatures (i.e., >300 • C) before~7 Ma. Although we do not have direct α-dose measurements for the exact same zircons that were dated with the ZHe method, we can use the spontaneous track density of the ZFT mounts of the exact same sample to obtain an approximate estimate of the likely mean α-dose for a given sample [45]. Future research would need the use of Raman spectroscopy on zircons (processed for ZHe dating) to obtain more direct information on accumulated α-dose. Nonetheless, the above interpretation is also supported by the fact that the analyzed zircons were colorless. The author of [46] showed that color in zircon is related to the accumulation of α-damage, but that color is lost when zircons are heated or reside at ambient temperatures of 325-475 • C, with full color resetting possible; therefore, annealing of α-damage is possible at temperatures as low as 350 • C [47]. Consequently, following the conclusions derived from the ZFT data, we can consider that the α-damage produced since crystallization has been annealed and only started accumulating since approximately <7 Ma. In that case, the oldest ZFT age at~7 Ma indicates that it should have begun to cool before this date. Cooling initiation at 7 Ma would induce an α-dose of approximately 6 × 10 15 to 4 × 10 16 α/g (Figure 7). Consequently, we propose that at a very low damage dose, the ZHe closure temperature is close to 80 • C (Figure 8). This result agrees with the prediction from [18] made using the theoretical approach of the diffusion behavior in zircon (Figure 8). This result has major implications because it demonstrates that the ZHe method has a large range of temperature sensitivities. A low closure temperature, similar to AHe, should be considered for zircons having a low damage dose (Figure 8). To further document and refine the models from field-based verification, it would be useful to integrate young plutons and/or young thermochronological data (few Ma; e.g., [48][49][50]) with as little inheritance as possible in such approach. sensitivities. A low closure temperature, similar to AHe, should be considered for zircons having a low damage dose (Figure 8). To further document and refine the models from field-based verification, it would be useful to integrate young plutons and/or young thermochronological data (few Ma; e.g., [48][49][50]) with as little inheritance as possible in such approach.  [44] to predict AHe, ZHe, AFT, and ZFT ages (in the blue box) using the model and model parameters of RDAMM [13], ZRDAAM [15,51] and the parallel curvilinear model of [51], respectively. This time-temperature path is an adaptation for higher closure temperature systems of the time temperature paths proposed for the southern Abancay Deflection by [27].  [44] to predict AHe, ZHe, AFT, and ZFT ages (in the blue box) using the model and model parameters of RDAMM [13], ZRDAAM [15,51] and the parallel curvilinear model of [51], respectively. This time-temperature path is an adaptation for higher closure temperature systems of the time temperature paths proposed for the southern Abancay Deflection by [27].
Minerals 2022, 12, x FOR PEER REVIEW 10 Figure 6. Plot of the Machu Picchu ZFT data (light red area) in the simplified spontaneous density ( s) and α-dose relationship from [45]. Grey dots are the 336 zircons they analyzed, a dashed and full lines represent the approximate and full relationships they computed, respec Spontaneous track densities from the Machu Picchu zircons indicate α-doses below 4 × 10 16 α Table 4. Closure temperature (Tc) estimation with a scenario implying rapid cooling initiat tween 5 and 7 Ma (case 1) and a simple scenario where the α-dose received is only due to the the Permo-Triasic, pluton 222 Ma; [23], (case 2) or of a >1 Ga pluton (case 3). Tc is estimated fo aliquot by reporting the computed α-dose for each case in Figure 8.       Red circles correspond to [15] data recalculated with the density functional theory results (see details in [18]), the yellow box places Pyrenean samples [52,53], and the blue box represents [10] highly damaged zircons. The black dashed line shows the shape of the closure-temperature/α-dose relationship presented in previous models. We added an estimation of the closure-temperature/α-dose relationship inferred from [18] (red dashed line), and the estimated closure temperature of the ZHe from Machu Picchu (this study, red square). Table 4. Closure temperature (Tc) estimation with a scenario implying rapid cooling initiated between 5 and 7 Ma (case 1) and a simple scenario where the α-dose received is only due to the age of the Permo-Triasic, pluton 222 Ma; [23], (case 2) or of a >1 Ga pluton (case 3). Tc is estimated for each aliquot by reporting the computed α-dose for each case in Figure 8.

Geological Implication
In a relatively old batholith such as the Machu Picchu Batholith, two geological processes could explain a low α-dose: (1) an important reheating >300 • C that reset all α-doses before~7 Ma could induce such a low α-dose, or (2) the batholith stayed at a relatively important temperature until recently, preventing any radiation damage effects until this date. Interestingly, previous studies [20,27], did not find evidence of any recent important reheating. Moreover, our ZFT data are not reset and show a low number of tracks. The author of [54] demonstrated that annealing of ZFT is much more temperatureand time-sensitive than healing of radiation damage in zircons. In our setting, zircon samples accumulated very little radiation damage in their lattice before a minimum of 7 Ma, suggesting that no important reheating occurred recently. We speculate that the batholith was at temperatures of~400 • C at 40 Ma ( Figure 5). We are not able to validate or disclaim if the Machu Picchu batholith experienced a reheating prior to 40 Ma with the dataset. Given that the intrusion age is~222 Ma, it is likely that the rocks experienced reheating at some point, but clearly before the rapid cooling starting at about~7 Ma. In such context, partial annealing of α-damage and full annealing of spontaneous fission tracks is very likely at temperatures of 300-400 • C and above e.g., [46]. Full annealing of α-damage is likely at temperatures above 530 • C. Without other higher temperature thermochronological data, we cannot well constrain that part of the thermal history. Nonetheless, our combined low-temperature thermochronology dataset does not indicate any reheating since~7 Ma. The simplest explanation, close to the time-temperature path used in Figure 5, is that the batholith stayed at a temperature higher than~300 • C until recently (before~7 Ma, oldest age of our dataset) before cooling below that temperature. AHe and AFT data in the southern part of the Abancay Deflection indicate a cooling acceleration at 5 ± 2 Ma from 100-150 • C [27]. Our ZHe and ZFT data show that this cooling acceleration may have initiated at least sometimes before 7 Ma and from a higher temperature (300 • C) than previously proposed. This suggests a higher cooling rate (>43 • C/Myr) than the one deduced from AHe and AFT data (21±6 • C/Myr) [20]. The AHe-and AFT-derived exhumation rates range between 0.6 and 1.9 km/Myr [20]. Taking into account ZFT data, it leads to an exhumation rate within the range of 1.3 and 3.0 km/Myr since~7 Ma, assuming end-member values for geothermal gradients of 26 ± 8 • C/km [55] or 18 ± 4 • C/km [27]. Our new data complete and strengthen the previous interpretations of the exhumation rates and further validate the tectonic syntaxis implications for the Abancay Deflection, as proposed in [27].

Conclusions
We provide new zircon fission track and zircon (U-Th)/He data for the Machu Picchu (Abancay Deflection, Peru) age-elevation profile. The new zircon ages are young (<7 Ma), reinforcing our previous interpretation of a young exhumation pulse (~5 Ma) and further favoring the interpretation of the Abancay Deflection as a tectonic syntaxis.
The ZHe system shows closure temperatures lower than those for apatite fissiontrack and apatite (U-Th)/He systems. This apparent contradiction, explained by low radiation damage in the zircons, is due to a simple time-temperature history with a longer stay at elevated temperature before cooling than to the age of the Machu Picchu batholith itself. In the present geological setting, these new data (1) demonstrate that the ZHe system is a thermochronometer not only sensitive to temperatures higher than 150 • C, but, in some cases, to temperatures lower than 100-110 • C or even 80 • C, and (2) complete the description of a setting (young ZHe ages, very low α-dose) that was predicted, but not observed by recent modelling. Our work highlights the importance of being aware of both the α-dose and the time since when damage accumulates to avoid biases when interpreting thermochronological data. In summary, our ZHe data open a door to designing further experiments and gaining a better understanding of He behavior in low radiation-damaged zircons.