Interpreting and reporting fission-track chronological data

The application of fission-track data to address a broad range of geological topics is based on the unique characteristic of fission tracks, which are produced by the spontaneous fission of ²³⁸U, to heal (i.e., fade naturally or thermally anneal) at elevated temperatures. Since each track is formed at a different time during the life of the host mineral, the fission-track method is a unique thermochronological tool in that each daughter product, which can be individually viewed microscopically, becomes a sensitive recorder of the thermal history of the host mineral. As such, daughter products frequently bear little or no relationship to the original formation age of the rock being studied. With the analysis of fission tracks in their host minerals, this process of thermal recovery—in tandem with the integration of relevant geological constraints—led to the development of robust modeling techniques for reconstructing the thermal histories of rocks and the rates of some geological processes in the upper crustal environment. Over the past 30 years, fission-track studies evolved into an important tool for thermochronology, providing opportunities for research over a broad spectrum of geological time covering a wide range of topics related to thermo-tectonic histories in different geological settings. Although thermochronological studies dominate, fission-track data are also occasionally applied as a geochronometer for dating discrete geological events. The basic principles and applications of fission-track dating were reviewed in several works (e.g., Fleischer et al., 1975; Wagner and Van den haute, 1992; Gallagher et al., 1998; Dumitru, 2000; Gleadow et al., 2002; Donelick et al., 2005; Tagami and O'Sullivan, 2005; Gleadow and Seiler, 2014; Malusà and Fitzgerald, 2019). For a historical perspective on the development of fission-track thermochronology, see Hurford (2019).

An early attempt at standardizing the treatment of fission-track data was proposed by Naeser et al. (1979), and further recommendations were made by Hurford (1990) to the International Union of Geological Sciences (IUGS) Subcommission on Geochronology, which also brought together different views on system calibration. Since 1990, major advancements in the methodology (e.g., Van den haute and De Corte, 2013) have necessitated the development of new, updated community agreed-upon data reporting schemas. These advancements include, but are not limited to, the measurement of an ever-growing suite of mineral compositional data and kinetic indicators (e.g., Carlson et al., 1999; Ketcham et al., 1999; Barbarand et al., 2003a), the increasing use of laser ablation–inductively coupled plasma–mass spectrometry (LA-ICP-MS) for ²³⁸U-content determinations (e.g., Hasebe et al., 2004; Chew and Donelick, 2012; Cogné et al., 2020) and elemental mapping (Ansberque et al., 2021), routine numerical thermal history modeling of fission-track data (e.g., Corrigan, 1991; Willett, 1997; Gallagher, 1995, 2012; Ketcham et al., 2000; Ketcham, 2005), deconvolution of detrital grain populations (e.g., Brandon, 2002; Dunkl and Székely, 2003; Jasra et al., 2006; Gallagher et al., 2009; Vermeesch, 2009, 2012, 2018), experimental characterization of very short-term annealing for fault-motion studies (e.g., Tagami, 2005, 2019), integration of fission-track data with other geochronological techniques on single crystals (double- or triple-dating; e.g., Danišík, 2019), and data visualization based on the production of large data sets (e.g., Kohn et al., 2002, 2005; Gleadow et al., 2002; Boone at al., 2021), as well as the advent of digital fission-track analysis (Gleadow et al., 2009). Consequently, large discrepancies remain in the detail and reporting of fission-track data and their associated analytical procedures, geosample metadata, geochemical analyses, and thermal history modeling protocols and results. These inconsistencies hinder public data transparency, accessibility, and reuse, and often impede Big Data regional syntheses and interlaboratory analytical comparisons from being readily performed.

The need for greater fission-track data transparency and community agreed-upon data reporting schemas are highlighted by other recent work, such as in the burgeoning field of machine learning-powered digital fission-track analysis (Nachtergaele and De Grave, 2021; Li et al., 2022), which will require ever-larger fission-track image data sets with systematic metadata descriptions to train artificial neural networks. Similarly, the development of the first bespoke relational data platform capable of storing, geospatially displaying, and interrogating fission-track analyses on a global scale (AusGeochem; Boone et al., 2022, 2023) requires analyses to be consistently reported in structured formats. In light of the aforementioned developments, and in response to roundtable discussions at recent international conferences on thermochronology, it is now timely to revisit and reset guidelines for reporting fission-track data. In this contribution, besides providing a brief background on the fission-track chronometer, methodologies, applications, and interpretations of data in different geological settings, we recommend guidelines for fission-track data reporting to be adopted by the entire thermochronological community, with the purpose of ensuring that they increasingly conform with FAIR (Findable, Accessible, Interoperable, and Reusable) data principles (e.g., Wilkinson et al., 2016; Stall et al., 2019; Devaraju et al., 2021; Klöcking et al., 2023) so that research data are increasingly useful, standardized, and machine-readable.

Fission tracks were first observed under a transmission electron microscope (TEM), appearing as linear cylindrical regions ~10 nm in width in mica (Silk and Barnes, 1959). Price and Walker (1962) demonstrated that such radiation damage trails in mica are highly reactive regions that could be “developed” and fixed permanently by etching them in hydrofluoric (HF) acid, which makes them readily visible at the micron-scale under an optical microscope and thus avoids track fading, a phenomenon commonly observed under the TEM. Price and Walker (1963) suggested that the spontaneous fission decay of ²³⁸U could form the basis for a new geological dating method for micas, and the suitability of other minerals and glasses was subsequently explored (Fleischer and Price, 1964a, 1964b). The simplicity of revealing fission tracks led to a rapid expansion of studies, which showed that they could be observed in a variety of natural minerals and glasses using different etching recipes; however, most of these minerals contain insufficient U or are too rare to be useful in routine dating studies (see Fleischer et al., 1975; Wagner and Van den haute, 1992). Pioneering studies by Naeser (1967, 1969) and Wagner (1968) focused attention on commonly occurring, U-bearing (typically in the parts per million to hundreds of parts per million range) accessory minerals such as apatite, zircon, and titanite, which proved to be particularly useful for geological studies. Apatite and zircon are by far the most commonly analyzed minerals in routine fission-track studies. The fundamentals laid down by these and other studies withstood the test of time, but numerous subsequent refinements define modern practice (e.g., Malusà and Fitzgerald, 2019).

For approximately every two million atoms of ²³⁸U that decay to ²⁰⁶Pb via a series of α and β reactions, one undergoes spontaneous fission. The nucleus splits into two, or occasionally three, smaller nuclei that are propelled away from each other. Each fission decay releases ~200 MeV of energy, of which ~169 MeV is kinetic energy imparted to the daughter nuclei, and the rest is released via neutrons, gamma rays, and subsequent decay of daughter products. As the daughter nuclei pass through the enclosing mineral, they lose their energy via ionic, thermal, and ultimately kinetic interactions, causing intense localized damage to the crystal lattice. The damaged regions, or fission tracks, are up to 9 nm in diameter and ~23 µm long in apatite, and somewhat shorter or longer in other minerals based on their density and atomic constituents (Paul and Fitzgerald, 1992; Jonckheere, 2003; Li et al., 2010). Note that spontaneous fission also occurs in ²³⁴U, ²³⁵U, and ²³²Th, but their spontaneous fission half-lives are too long and/or abundances too low to produce a significant number of natural tracks compared to those produced by ²³⁸U.

As indicated above, to be observed directly using optical microscopy, fission tracks must be etched (Price and Walker, 1962). The etching velocity along the track (v_T) is much faster than through the bulk grain (v_B). Thus, an etched track has a diagnostic linear appearance created by dissolution of the ultrathin damage zone, followed by slower etching into the surrounding, more pristine lattice. Some of the latent track is not revealed by etching, as etching rates are not enhanced toward the initially formed tips, which results in a length deficit (Jonckheere, 2003). This is in no way debilitating to the method, but it emphasizes that virtually all theory and practice of fission-track thermochronology is in reference to the etchable portion of damage.

Sample preparation consists of embedding mineral grains in epoxy or Teflon (depending on the etchant required) and polishing the mount to expose an internal surface at least 20 µm into the grains so that the center of the exposed grain has no exterior-intersecting tracks present (Figs. 1A and 1C). After the mount is etched, two principal measurements are made on crystals oriented with the crystallographic c-axis (for apatite and zircon) parallel to the polished surface; these are made under transmitted and/or reflected light (Figs. 1B and 1C). Ages are determined from track density (N/cm²), which is the number of tracks that intersect the polished grain surface per unit area (Fig. 1C). Confined track lengths are revealed when etchant travels into the grain interior via surface-intersecting pathways and intersects tracks at depth that have both tips within the solid grain. Etch figures, the pits formed by a track descending into the polished apatite grain surface, are also commonly measured; they are elongated parallel to the c-axis, and their corresponding diameter, D_par (Fig. 1D), can be used as a kinetic indicator for both initial (unannealed) track length and relative annealing resistance in apatite (Burtner et al., 1994).

The potential time span covered by fission-track dating covers a broad range. Ages are reported for man-made objects (doped with U) dating to several tens of years, mineral inclusions in heated rocks, baked soils, pottery dating to several thousands of years, as well as minerals in rocks and meteorites older than 1 Ga (Fleischer et al., 1975; Wagner and Van den haute, 1992). The age range covered is principally governed by a sample's areal density of fission tracks (tracks/cm²) and assumes that materials of interest have adequate age or adequate U content, or both. To calculate a precise age, if possible, hundreds—or where feasible, thousands—of tracks should be counted. However, this is not always possible, particularly in young, low-U grains in which useful data may still be acquired even though uncertainties will be large. The age limits imposed by fission-track dating are related to the product of age and U content and the somewhat subjective consideration of the time taken to conveniently count enough tracks for a robust age determination. Figure 2 shows the typical dateable age limits using optical microscopy for the most commonly dated minerals and glasses in relation to their U content and fission-track density.

A fission-track age represents the integrated record of a mineral grain's thermal history starting from the time it last passed into a thermal regime in which fission tracks are partially or fully retained. Therefore, a simplistic interpretation of an age as corresponding to an original formation age (or event), such as the time when the sample cooled to below the system closure temperature, is only true in a limited number of cases, such as volcanic eruptions (e.g., Westgate, 1989; Gleadow et al., 2015), kimberlite and diatreme emplacement (e.g., Brookins and Naeser, 1971), localized conductive heat transfer related to igneous activity (e.g., Calk and Naeser, 1973; Tagami and Shimada, 1996), shock-wave heating during meteorite impact (e.g., Storzer and Wagner, 1977), and archaeological studies (e.g., Wagner, 1978; Kohn, 2017). Instead, data are mostly interpreted in terms of the ways they elucidate thermal histories, often with the aid of computational tools. As such, fission-track data can be used to quantify the timing, rates, and distribution of a breadth of geological processes that can affect the thermal state of the crust over geological time. These include the advection of mass and heat due to the growth of mountain belts, extensional basin and passive margin formation, fault zones, and long-term landscape exhumation (e.g., Ehlers, 2005; Stockli, 2005; Reiners and Brandon, 2006; Fitzgerald and Malusà, 2019; Schildgen and van der Beek, 2019; Tagami, 2019; Wildman et al., 2019; Kohn and Gleadow, 2019). As such, fission-track thermochronology has been used to study a range of natural phenomenon through deep-time, including paleoclimate (e.g., Kohn et al., 1992; Shane et al., 1995), surface weathering processes (e.g., Ault et al., 2019), and climatic change (e.g., Herman et al., 2013), as well as to constrain the formation and preservation of various natural resources, such as hydrocarbons (e.g., Schneider and Issler, 2019), and hydrothermal and supergene ore deposits (e.g., McInnes et al., 2005; Gong et al., 2021; Sun et al., 2021). In certain instances, fission-track analyses can also record groundwater advection (e.g., Boone et al., 2016), hydrothermal fluid flow (e.g., Duddy et al., 1994), and wildfires (e.g., Reiners et al., 2007).

Studies involving the use of different etchants revealed the presence of fission tracks in more than 90 U-bearing minerals and glasses (Fleischer et al., 1975; Wagner and Van den haute, 1992). However, their abundance, U content, and track-stability traits result in very few being routinely analyzed for fission-track thermochronology. The most widely dated minerals are apatite, zircon, and titanite; but more recently, the characteristics of monazite suitability for fission-track studies were further evaluated (e.g., Jones et al., 2021). These are all commonly present to different degrees as accessory minerals in many igneous and metamorphic rocks and as detrital components in some sedimentary rocks. Volcanic, impact, and man-made glasses and pseudotachylites were also occasionally studied. For a summary of preferred lithologies hosting target minerals for fission-track studies and a general guide for field collection, see Kohn et al. (2019).

Fission-track ages are calculated using the equation for branching decay age:

where P and D are the respective numbers of parent atoms and daughter products (i.e., fission tracks) as measured today after some time, t, has elapsed; λ and λ_f are the total ²³⁸U-decay constant (1.55125 × 10⁻¹⁰ yr⁻¹; Jaffey et al., 1971) and the fission-decay constant (7.9–8.7 × 10⁻¹⁷ yr⁻¹; Holden and Hoffman, 2000), respectively. Spontaneous fission decay of ²³⁸U is measured as the density of tracks crossing an internal counting surface (ρ_s). The latent track density per unit volume (D) is estimated from this as:

where R_s is the average etchable range of a single spontaneous ²³⁸U fission-track fragment (in μm), and ηq reflects detection efficiency, which is a combination of physical (η) and analyst-specific (q) factors that are difficult to disentangle (Fleischer et al., 1975; Jonckheere and Van den haute, 1996).

For fission-track dating, the most common strategy traditionally employed for studying minerals has been the external detector method (EDM), which involves sending off polished grain mounts for thermal neutron irradiation in a nuclear reactor. However, more recently, laser ablation–inductively coupled plasma–mass spectrometry (LA-ICP-MS) emerged as an alternative method (Svojtka and Košler, 2002; Hasebe et al., 2004). An outline of the experimental procedures for these two methods and a comparison of their pros and cons is presented in Section 4.

For the EDM, an external fission particle detector (usually low-U muscovite) is attached to the polished mount surface, and the sample package is irradiated to induce fission of ²³⁵U, sending newly formed tracks into the detector, which can then be etched and counted. The resulting parent equation is:

where ρ_i is the induced track density; σ is the cross section of ²³⁵U for fission induced by thermal neutrons; ϕ is the neutron fluence received in the reactor; I is the natural ²³⁵U/²³⁸U ratio (1/137.8); ηq_det is the detection efficiency of the analyst for the detector; R₂₃₅ is the etchable range of a single ²³⁵U fission fragment (in μm); and g is a geometry factor (0.5) to reflect that tracks only originate from outside the detector, or with a 2π geometry, as opposed to the 4π geometry of the mineral internal surface where daughters are measured.

These daughter and parent terms include factors that can be difficult to know or measure with confidence, including counting efficiencies for grains and detectors, reaction cross sections that vary with reactor conditions, neutron fluxes, and absolute fractionation factors. Even the ²³⁸U fission-decay constant has historical values that vary by ~20% (Hurford and Green, 1981), although a value of 8.46 ± 0.06 × 10⁻¹⁷ yr⁻¹ (Spadavecchia and Hahn, 1967; Bigazzi, 1981) was used for age determinations using an absolute LA-ICP-MS approach (e.g., Hasebe et al., 2004). The practical solution was to use a zeta factor to calibrate all of the problematic variables against a set of age standards (Hurford and Green, 1983; Hasebe et al., 2004; Cogné et al., 2020).

For EDM dating, the neutron flux can be inferred by co-irradiating a U-doped standard glass with known, uniform composition. By measuring fission tracks generated within that glass that are implanted in an external detector, one can estimate fluence as ϕ = Bρ_d, where B is a proportionality constant and ρ_d is the fission-track density in the detector. ρ_d = ρ_i/A_d is the number of induced fission tracks (ρ_i) counted over an area, A_d, in a second external detector that is attached to a co-irradiated dosimeter glass. Zeta (ζ) can then be defined as:

as given by Hurford and Green (1983), although this definition assumes that the length and efficiency terms cancel, an early assumption (Fleischer et al., 1975) that is not strictly correct (Jonckheere and Van den haute, 1998), and so a full consideration of zeta should include them as well.

The resulting simplified age equation for the EDM, where t is the fission-track age (in Ma), is:

ζ can be calibrated by performing a series of measurements on age standards, whereby:

For LA-ICP-MS analysis, the analysis of an unknown is derived from an appropriate isotopic ratio (e.g., ²³⁸U/⁴³Ca or ⁴⁴Ca for apatite, ²³⁸U/²⁹Si for zircon) that is calibrated and corrected for fractionation and instrumental drift by analyses of an appropriate reference material with known ²³⁸U, ⁴³Ca, ⁴⁴Ca, or ²⁹Si concentrations. An example for apatite is:

where ²³⁸U_r and ⁴³Ca_r are the known ²³⁸U and ⁴³Ca concentrations in parts per million of a reference material, ⁴³Ca_s is the known ⁴³Ca concentration of apatite assuming stoichiometry, m is the measured signal ratio (²³⁸U/⁴³Ca), and subscripts s and r refer to the sample and reference material, respectively. κ denotes a factor to correct the difference in elemental behavior during LA-ICP-MS measurement, and mainly depends on differences in ablation efficiency between sample and reference material, and fractionation between U and Ca through ablation, transport, and ionization in the plasma. Conventionally, the ²³⁸U concentration of the reference material is given in parts per million (μg/g); therefore, the resultant ²³⁸U concentration of the mineral is obtained in parts per million. Then, the number of parent atoms per unit volume is calculated as:

where N_A is Avogadro's number, d is the specific density of the dated mineral to convert mass to volume, and M²³⁸ represents the atomic weight of ²³⁸U, yielding the age equation for absolute LA-ICP-MS dating:

where N_s is the number of spontaneous tracks counted over an area A_s (in μm²), which is equivalent to ρ_s, g is a geometry factor (that = 1 in the case of LA-ICP-MS fission-track analysis if N_s was counted on an internal crystal surface), [²³⁸U] is the number of ²³⁸U atoms per μm³ measured by LA-ICP-MS, and ηq_s is a track-detection efficiency factor that depends on the mineralogical, etching, and observation conditions. (U can also be measured by other techniques such as secondary ion mass spectrometry (SIMS) or electron microprobe analysis, but LA-ICP-MS is assumed in the remainder of this section).

It is necessary to measure or reasonably assume values in these equations. The Ca or Si concentration and specific density of a mineral may vary among samples, and additional labor would be required to measure these values, which would also introduce additional sources of error. To overcome this problem and the others outlined above, the zeta approach is commonly used for determination of the parent ²³⁸U by the LA-ICP-MS dating method with the following general age equation:

where ζ is a calibration constant (in Myr μm²) obtained by analyzing a reference material of known t, and [U] is the U concentration (in parts per million) or the U/Ca ratio (denoted as P in Section 4.3). But the components constituting the ζ value in this case differ from those incorporated in the EDM approach in Equation 5, such that ξ (Xi) was proposed as an alternative aggregate constant to explicitly differentiate among them (Gleadow and Seiler, 2014; Seiler et al., 2023). Also note that for the LA-ICP-MS methods, either absolute dating or zeta-type calibration approaches are both possible (e.g., Hasebe et al., 2004, 2013; Chew and Donelick, 2012; Soares et al., 2014; Cogné et al., 2020; Seiler et al., 2023; see Section 4.3).

The spontaneous fission of ²³⁸U is a relatively rare event, and track counts are typically in the ones to tens, resulting in relatively large single-grain age uncertainties, especially in samples with low U content and/or relatively young cooling ages. Hence, for most studies it is necessary to combine observations from many grains in a single sample to overcome this imprecision. The likelihood that all grains analyzed in a sample comprise a single population, their analytical precision, and degree of age dispersion can be quantified statistically and visualized (Galbraith, 2005; Vermeesch, 2019). Therefore, the way in which observations are combined to determine an age and uncertainty can be approached in different ways. See Sections 3.2 and 3.3 for further discussion and Section 5 for possible causes of excess age dispersion.

The uncertainty budget of the fission-track method includes both random and systematic components:

1. Random errors affect all grains independently. They include the Poisson uncertainties of the fission-track counts and the effect of mass spectrometer detector noise on the precision of the U measurements.

2. Systematic errors affect all grains simultaneously, resulting in correlated uncertainties. They include the uncertainty of the ζ-calibration factor, the dosimeter track density ρ_d, the decay constants λ and λ_f, and the reproducibility of the internal standards used to compute the U concentrations or U/Ca ratios.

The standard errors of fission-track ages (s[t]) are obtained by standard error propagation:

For EDM (Equation 5):

For LA-ICP-MS absolute dating (Equation 9):

where the systematic uncertainties are marked by curly brackets, and it is assumed that the uncertainty of the U concentrations only accounts for the detector noise, and not for any internal or external calibration factors.

For LA-ICP-MS ζ-calibration dating (Equation 10):

The conventional way to deal with random and systematic uncertainties is to separate them in a so-called hierarchical error propagation approach. This means that the error propagation proceeds in two steps. The single-grain age errors are calculated first, using just the random sources of uncertainty. The systematic errors are only added later, for example, when the single-grain ages are averaged or subjected to some form of deconvolution algorithm (see Section 3.3). The systematic sources of uncertainty that should be included in the uncertainty budget depends on the purpose of the results. For example:

Compare a single-grain age estimate with an independent stratigraphic age constraint.

Strategy. Include all sources of uncertainty in Equation 11 or 13.

Compare the mean age of multiple grains with an independent age constraint.

Strategy. Assuming EDM data, (1) compute the single-grain age errors using only the random sources of uncertainty (i.e., N_s and N_i); (2) average the data, for example, using the random effects model of Section 2, and estimate the standard error of that average; (3) augment the standard error of the average with the systematic uncertainties associated with the dosimeter track density and ζ-calibration constant:

where s[t] and s[t] are the standard errors of the average age (t) with and without systematic uncertainties, respectively.

Assess the age difference between two samples that were jointly irradiated and counted by the same analyst.

Strategy. Compute the single-grain ages and average the data as in Scenario 2, but do not add the systematic uncertainties afterward.

The fission-track uncertainty budget is dominated by the spontaneous track count, N_s, which follows a Poisson distribution. Poisson distributions are completely described by a single parameter, λ_s, which is both their mean and their variance. λ_s is unknown but can be estimated from the data using the so-called method of maximum likelihood (Galbraith, 2005). For a single grain estimate, in which only one measurement is available, the maximum likelihood estimate of λ_s is simply N_s. Given the fact that the variance of a Poisson distribution equals its mean, the standard error of N_s is . Consider the sequence of hypothetical (but realistic) spontaneous track counts and their standard errors presented in Table 1. The uncertainties quickly grow with decreasing N_s. Consequently, the single-grain age uncertainties of fission-track data are much larger and more variable (“heteroscedastic”) than those of most other geochronometric methods. The only way to improve the precision of the method is to analyze multiple grains from the same sample and somehow average them. To assess the best strategy to do so, it is important to first subject the single grain analyses to a visual inspection. The radial plot is a graphical device that was specifically designed for this purpose (Galbraith, 1990); see Figure 3.

Given a set of estimates (z) and their standard errors (s[z]), the radial plot (Fig. 3) is a bivariate scatter plot that sets out against , where z_∘ is a reference value (e.g., an average value). Thus, the z-value of each aliquot scales with the slope of a line that connects it to the origin, with positive and negative slopes marking values that are greater and smaller than z_∘, respectively. The precision of each aliquot is given by its horizontal position on the radial plot, so that precise measurements plot to the right of imprecise measurements.

For nonnegative data such as single-grain age estimates and ratios, it is customary for the z-values to represent transformed versions of the data rather than their raw values. Examples are the logarithmic (for strictly positive values), square root (in the presence of zeros), and arcsine (for count ratios) transformations. Their main purpose is to reduce or eliminate the inherent skewness of nonnegative data sets, and thereby allow visual inspection over the entire range of values.

Projecting a “2-sigma” error bar around any given aliquot onto an appropriately transformed radial scale yields an approximate, asymmetric 95% confidence interval for the corresponding value (gray wedge in Fig. 3A). The radial plot also provides an easy way to visualize the degree to which the analytical uncertainties are consistent with the observed scatter of the data. If analytical uncertainty is the only source of scatter, then ~95% of the data should plot within a symmetrical 2-sigma confidence band around the origin (dashed lines in Figs. 3A–3D). Data sets that scatter significantly beyond this region are “over dispersed” with respect to the known analytical uncertainties and reveal either some additional measurement error or are indicative of geological complexity, which may be manifested as continuous or discrete age mixtures.

The degree of overdispersion can be formally assessed by defining a chi-square statistic (χ²). For the EDM, this number can be calculated directly from the raw data:

with and , where N_sj and N_ij are the spontaneous and induced fission-track count of the j^th grain, respectively (for 1 ≤ j ≤ n).

For LA-ICP-MS–based fission-track geochronology, it is not possible to quantify the dispersion from the raw data, but an equivalent statistic, (), can be defined using the log-transformed ages and their standard errors (Galbraith, 2010; Vermeesch, 2017):

The expected distributions of the chi-square statistic scales with (1) the degree of dispersion relative to the analytical uncertainties and (2) the size of the data set (n). To reduce the sample size dependency of χ² (or ), it is useful to divide it by the number of degrees of freedom (i.e., df = n − 1). The resulting “reduced chi-square statistic” is more widely known as the “mean square of weighted deviates” (MSWD) in geochronology. It is widely used in U-Pb and Ar-Ar geochronology but historically less so in fission-track analysis. Nevertheless, Figure 3 uses the MSWD to explore several dispersion scenarios.

Figure 3A shows synthetic EDM data that scatter beyond the 2-sigma confidence band of the radial plot. Under the “null hypothesis” that the analytical uncertainties are the only source of scatter, the probability (“p-value”) of observing the corresponding χ²-value of 154 (or a value greater than that) is only 1 × 10⁻²². This is extremely unlikely, resulting in a rejection of the null hypothesis, and leading to the conclusion that the data are significantly overdispersed.

Figure 3B shows a second data set whose measurements all fall well within the 2-sigma band. This high degree of clustering is unlikely to happen by chance. The probability of observing an χ²-value of 2.85 (or less) is only 1 × 10⁻⁵ (corresponding to a p-value of 0.999993). This data set is underdispersed with respect to the analytical uncertainties.

Figure 3C shows a data set for which all but one aliquot fits within the 2-sigma band. In the absence of overdispersion or underdispersion, the probability of observing an χ²-value that exceeds the observed value of 21 is 0.34. This p-value is neither extremely high nor low and provides no strong evidence for either overdispersion or underdispersion. Data sets that are neither underdispersed nor overdispersed are easy to interpret. They indicate that the sample cooled through a relatively small window of time that can be pinned down (for EDM data) using the pooled age, i.e., by treating N_s• and N_i• as if they were counted in a single grain.

Figure 3D shows a data set that was drawn from exactly the same distribution as the example shown for Scenario 3. The only difference is the order of magnitude increase in sample size (n = 200 versus n = 20). The degree of scatter in this radial plot is not noticeably larger than in Figure 3C, and the MSWD is of a very similar magnitude. However, despite this similarity in observable scatter, the corresponding p-value is two orders of magnitude smaller. Thus, even though the amount of overdispersion is minor, its statistical significance is high.

Underdispersion may be caused by incorrectly propagated analytical uncertainties, which may be caused by observer bias, or they may be caused by the inclusion of systematic effects in the error propagation. Alternatively, underdispersion may also indicate the violation of the parametric assumptions behind the chi-square test. This is most pertinent for the LA-ICP-MS approach, whose error distribution may not fulfil the lognormal assumption.

Underdispersion is usually avoidable and always undesirable. Overdispersion is neither of these things. Figure 3D shows that overdispersion is unavoidable given a large enough sample size. Furthermore, the presence of overdispersion is neither good nor bad, although it does somewhat complicate data interpretation. In some cases, instead of registering a single, discrete geologic event, overdispersed data sets record more complex geologic histories. How to process the data then depends on the nature of those histories (Fig. 4).

Instead of assuming that all grains in a sample register exactly the same discrete cooling event, fission-track data can also be modeled assuming a continuous (lognormal) age distribution. The “random effects model” aims to constrain this distribution in the presence of (Poisson or lognormal) measurement error. The central age corresponds to the mean of the (log-transformed) continuous age distribution, and the “overdispersion” is given by its standard deviation. The random effects model is a useful way to capture the effects of compositional variability (e.g., variable Cl⁻ content in apatite) on the annealing temperature of slowly cooled samples. The central age and overdispersion parameter parsimoniously describe the diachroneity of the system under this scenario.

The random effects model is only appropriate when all grains in a sample experience the same thermal history. It is generally not applicable to detrital samples, unless they were completely annealed after deposition. Non-reset or partially reset sediment samples may contain mixtures of discrete age components that were sourced from different catchment areas. Deconvolution of these mixtures is possible under the assumption that the data scatter around them is solely due to analytical uncertainty.

It is not always easy to decide how many components should be fitted to the data. Automated model selection rules, such as the Bayes information criterion (BIC; Galbraith, 2005), should be used with caution as the number of components that they return always increases with sample size. The main purpose of automated mixture modeling is to provide a maximum limit on the number of statistically distinguishable components. In other words, if the BIC returns three components, then it is definitely not justified to constrain four.

Finite mixture modeling is usually not a good way to constrain the maximum depositional age of a sedimentary deposit. This is because the sample-size dependency of model complexity causes the youngest age component to drift to younger values (see fig. 6.3 of Vermeesch, 2019, for an example). The minimum age model of van der Touw et al. (1997) is a more appropriate solution that combines aspects of the random effects and finite mixture models. It assumes that the observed age distribution combines a discrete youngest age component with an older lognormal distribution that is truncated at said youngest age value. In contrast with the BIC and other automated model selection approaches, the minimum age model does converge to a distinct minimum age.

In principle, mixtures of continuous distributions are the most flexible way to model complex age distributions (Jasra et al., 2006). However, in practice, such models are rarely useful. For a k-component system, mixtures of random effects models involve 3k − 1 parameters (k central ages, k overdispersion parameters, and k − 1 proportions). Searching such a high-dimensional parameter space is not only computationally demanding and unstable, but also requires extraordinarily large data sets for any k > 2, for example.

Following field collection and conventional mineral separation, mineral separates are mounted in epoxy or Teflon, followed by grinding and polishing to expose flat, internal mineral grain surfaces. Samples are then etched to reveal tracks using protocols particular to the mineral studied (e.g., Figs. 1B and 1C). For further details on sample preparation and chemical etching schemes for fission-track analysis of commonly dated minerals and glasses, see Kohn et al. (2019); for monazite, see Jones et al. (2019).

Several strategies for acquiring fission-track ages have been described, each consisting of a specific sequence of steps (e.g., Hurford, 2019). As indicated earlier, only the mainstream EDM and LA-ICP-MS methods will be considered here. A first step common to both of these methods is to examine the etched grains using a research-grade microscope and select grains suitable for further examination. Criteria for optimal microscope viewing, grain selection, and recording of grain position are outlined in Kohn et al. (2019). Once this information is recorded, each selected grain is examined, and the number of spontaneous fission tracks is counted. Sample track density is calculated, as required for the age equations (see Section 3.1). Spontaneous fission-track counting and the acquisition of other important fission-track parameters, such as confined track-length and D_par measurements, can be accomplished by two different approaches. The conventional method, in which measurements are acquired manually, can be labor intensive. In recent years, an alternative procedure utilizing autonomous digital microscopy for image capture and subsequent digital image analysis for the recognition and counting of fission tracks in minerals (e.g., apatite) was developed (e.g., Gleadow et al., 2009, 2019). For further discussion here, we differentiate between these two strategies by using the terms “manual” and “image-based.” Note that for the latter approach, the application of a thin, metallic coating, typically gold, is routinely applied to the etched grain surface prior to image-based analysis. This eliminates strong internal reflections during image capture resulting from the low reflectivity of polished grains (particularly apatite), which diminishes the quality of the complementary reflected light image required for the automatic “coincidence mapping” fission-track identification routine (Gleadow et al., 2019). The workflows, in terms of steps employing these two approaches for apatite fission-track (AFT) dating of a batch of 10–12 samples, are compared in Figure 5, with the approximate time involved. Further details are outlined below in various sections.

The sequence of steps involved in preparing samples of different minerals for irradiation, and details of measurement strategies using the EDM, were described in several works (e.g., Gleadow et al., 2002; Tagami and O'Sullivan, 2005; Kohn et al., 2019). Briefly, spontaneous fission tracks are counted on an internal polished surface of etched mineral grains and induced fission tracks in an external track detector (usually a low-U muscovite) firmly in contact with the mineral surface. Following irradiation, the mica external detectors are etched to reveal induced tracks corresponding to each grain in the mineral mount. As the grain mount is not re-etched postirradiation, only spontaneous tracks are observed in sample grains, while only induced tracks are revealed in the external detector. The etched surfaces therefore have a mirror-image relationship to one another and are mounted side-by-side on a microscope slide so that precisely matched areas of spontaneous and induced tracks can be counted. Using coordination markers, the matching between grains and mica is usually automated by using a computer-operated stage. An EDM age determination therefore requires selecting suitable grains in the mount, counting the spontaneous tracks present, and then counting induced tracks over the corresponding area on the external detector. Poorly oriented grains (e.g., apatite, zircon, or titanite grains whose polished internal faces are not parallel to the c-axis) and those with nonuniform U distributions, dislocations, inadequate polishing, or other defects are excluded.

During neutron irradiation in a well-thermalized reactor, it is essential to monitor the total neutron fluence received by each grain mount to establish how much ²³⁵U underwent fission during irradiation. This is achieved by including neutron flux monitors in the same package as the sample mounts. These are also usually muscovite external detectors attached to reference glasses with a known and homogeneous U distribution, which allows glass-induced track densities to be determined in the etched mica detector. The monitors are placed at the top and bottom of the irradiation package so that the total neutron dose received can be monitored and any flux gradients corrected for. Once the neutron fluence is known for each sample in the irradiation package, U concentrations can be calculated for each mineral grain analyzed based on the corresponding number of induced fission tracks recorded in the external detector. Although flux monitors may be calibrated against the absolute thermal neutron fluence, the fluence calibration is usually included in the zeta constant so that the standard glass track density is applied directly in the age calculation (Equation 5). With the progression of various methodologies, however, some research groups pursued absolute EDM dating, using measurements of efficiency and flux and agreed-upon values for the cross section and decay constants (Jonckheere, 2003; Enkelmann et al., 2005).

Most LA-ICP-MS analytical protocols (Hasebe et al., 2004, 2013; Chew and Donelick, 2012; Soares et al., 2014; Cogné et al., 2020, Cogné and Gallagher, 2021; Seiler et al., 2023) undertake fission-track dating by spot ablation, rather than the more time-consuming multispot analyses or mapping approaches that use rastering (line scans). LA-ICP-MS spot analysis also easily facilitates the simultaneous acquisition of a U-Pb age (Chew et al., 2014a), petrogenetically diagnostic trace elements (e.g., Sr, Y, and the rare earth elements (REEs) for protolith fingerprinting of apatite detritus; O'Sullivan et al., 2020), and elemental concentrations that are kinetic parameters in fission-track annealing models (e.g., analysis of Cl in apatite; Chew et al., 2014b). In the following sections, it is assumed that simultaneous U-Pb age information is also acquired as it has implications for the analytical protocol employed.

Instrument tuning should be performed before each LA-ICP-MS session on a well-characterized reference material, such as the National Institute of Standards and Technology (NIST) 612 doped glass, to minimize fractionation for U-Pb geochronology, to produce Th/U ratios close to the atomic ratios of the glass (~1.01), and to produce low oxide production rates (ThO⁺/Th⁺ < 0.15%) while optimizing the signal intensities of ²³⁸U and an element that is assumed to be stochiometric (⁴³Ca or ⁴⁴Ca in apatite or ²⁹Si in zircon). With the fast washout (i.e., low dispersion) of modern laser cells, between 5 s and 10 s of baseline acquisition suffices. The pulse energy should be set comfortably above (~20%) the ablation threshold on the LA-ICP-MS system for NIST glass and apatite or zircon, which minimizes evaporation of volatile elements during ablation, and should yield ablation depths of ~15 μm for 250 shots (analysis duration of 25 s for a laser repetition rate of 10 Hz). This ablation depth generates sufficient analyte volume for precise simultaneous acquisition of U-Pb isotopes and trace-element analyses. A laser ablation spot size of 30 μm or 35 μm is typically employed for LA-ICP-MS apatite analyses, which yields sufficient signal for U-Pb geochronology while not yielding too high a pit aspect ratio, which would exacerbate downhole fractionation in U-Pb geochronology. For zircon, the spot size may need to be smaller (~20 μm) if U-Pb geochronology is being undertaken to keep the ²³⁸U signal below the pulse-analog (P/A) switchover threshold (typically between 1 Mcps and 4 Mcps).

A standard reference material glass (NIST SRM 610 or 612) is typically used as the primary reference material to correct for LA-ICP-MS session drift. Normalization relative to a stoichiometric internal standard isotope in the reference material removes much of the effects of session drift; nevertheless, minor residual drift in Ca- or Si-normalized elemental abundances may still remain and necessitates drift correction. The ²³⁸U/^43,44Ca or ²³⁸U/²⁹Si ratio can be depth weighted using a spherical depth-weighting function (Chew and Donelick, 2012) so that U-concentration data close to the grain surface are weighted more heavily than U concentrations at depth, down to the distance of one fission fragment range (e.g., ~8 μm for apatite) from which any surface-intersecting track can arise. U zoning with depth can thus be easily accounted for using this approach, and the comparison of the depth-weighted and raw ²³⁸U/^43,44Ca or ²³⁸U/²⁹Si ratios from the same ablation is a useful method for identifying the presence of U zoning (Cogné et al., 2020). However, lateral U zoning on the etched surface of the grain mount (i.e., the counting area) is more difficult to assess, especially for low-track densities, and is discussed in Section 4.4.

As mentioned previously, absolute dating using the EDM approach was advocated by Jonckheere (2003) and Enkelmann et al. (2005). LA-ICP-MS fission-track dating can potentially also be used as an absolute dating technique (Equation 9), based only on best estimates for various constants in the age equation as an alternative to the more conventional zeta approach, which is dependent on calibration against independently dated age standards. Since the pioneering studies of Hasebe et al. (2004, 2009), little work has been reported to test an absolute dating approach in terms of calibration and standardization. The most detailed study since those works was reported by Seiler et al. (2023), who presented an extensive data set of 572 apatite fission-track single-grain ages from 17 samples that represent a wide range of fission-track ages, and U content and track densities, which were analyzed using both the EDM and LA-ICP-MS techniques. Both approaches were applied to the same selected grains (whose EDM micas were also screened for possible U zoning) using identical spontaneous track densities. A strong 1:1 correlation was observed among the single-grain, fission-track ages (within analytical uncertainties) analyzed by both methods, irrespective of sample age or cooling history (Fig. 6). However, inter-grain variability of ²³⁸U was significant, as the LA-ICP-MS determinations were more precise, which may explain the greater dispersion of single-grain ages calculated by that method compared to those obtained by the EDM.

Despite the good agreement observed in single-grain data reported by Seiler et al. (2023), pooled and central ages determined by both methods did not always correlate so well. The observed age difference between techniques appears to be negatively correlated to mean track length. EDM- and LA-ICP-MS–derived ages for rapidly cooled, well-established apatite reference materials such as Durango or Fish Canyon Tuff, with long mean track lengths (>14 µm), were statistically insignificant, and the corresponding pooled/central ages that correlate within error were statistically indistinguishable. However, in other samples with more complex cooling histories, a significant difference of up to 25% was observed among central ages determined by LA-ICP-MS and EDM, with the latter being consistently younger. By contrast, samples with moderate to shorter mean track lengths (<12–13µm) yield LA-ICP-MS pooled/central ages that are consistently older than the EDM ages, which reflects the degree of track annealing experienced by those samples with slower, more protracted or complex cooling histories. This discordance is attributed to the fact that fission-track ages derived by EDM use the almost universally common zeta-calibration approach, which implicitly assumes that the mean track length of a sample (2R_s) is equal to that of the rapidly cooled age standards from which it is derived. However, for the absolute LA-ICP-MS approach, the etchable range of a single spontaneous fission-track fragment, R_s, is a fundamental component (Equation 9), and the calculated ages vary with their mean track lengths, but such an equivalent term is not included in the EDM zeta age equation. Hence, any fission-track age that uses the zeta-calibration approach (for both EDM and absolute LA-ICP-MS–based data) implicitly assumes that the mean track length of the sample equals that of the age standard. Hence, the two approaches only produce concordant results when they are based on similar assumptions regarding the mean etchable range of fission fragments (i.e., rapidly cooled samples with long mean track lengths), and as such most fission-track ages calculated in this way should be regarded as model ages (Seiler et al., 2023).

Chew and Donelick (2012) and Cogné et al. (2020) employed a zeta-based approach to LA-ICP-MS apatite fission-track dating, using the equations presented in Donelick et al. (2005). While these studies used crushed Durango apatite (sieved to 200–300 μm shards) as the primary zeta reference material, it is also possible to use a large c-axis–parallel slab from a single Durango crystal for analysis of individual sub-areas. As Durango apatite is gem-quality and largely devoid of defects, counting using a 50× objective is possible, which maximizes the size of the area to be counted (A_s) and thus the number of spontaneous fission tracks (N_s) counted. Cogné et al. (2020) analyzed 80–100 shards of Durango during a large, primary LA-ICP-MS session, which was repeated two more times, and then the U/Ca ratios of each shard were averaged. Each counted shard is large enough (~2 × 10⁻⁴ cm²) to subsequently accommodate between 15 and 25 laser spots during the subsequent LA-ICP-MS session with unknowns. The pool of 80–100 Durango shards allowed for ~7000–10,000 separate single-grain, fission-track age determinations (assuming four Durango zeta shard analyses for every 20 single-grain unknowns). Once the Durango zeta material is used up, the mount can be gently repolished down by a depth of ~20 μm to remove the ablation pits of previous analyses and then be reused. While ²³⁸U/⁴³Ca ratios vary from one shard to another, it is assumed that the ²³⁸U/⁴³Ca ratio on the scale of the counted shard is constant based on (1) LA mapping of the Durango zeta standard prior to crushing, (2) the consistency between depth-weighted ²³⁸U/⁴³Ca ratios and non-depth–weighted ²³⁸U/⁴³Ca ratios, and (3) the internal homogeneity with mean standard deviations of ~2.0% (maximum of ~5%) among the three analyses of each shard in the primary zeta session (Cogné et al., 2020).

Cogné et al. (2020) employed a primary zeta LA-ICP-MS session to enable calculation of the primary zeta factor, ζ_ICP, and its error, s_ζICP. The primary zeta factor, ζ_ICP, can then be employed in all subsequent LA-ICP-MS sessions involving the analysis of unknowns. This avoids the time-consuming counting and analysis of many Durango shards for every secondary LA-ICP-MS session that would otherwise be required to yield a relatively precise zeta factor. A session-specific zeta fractionation factor (X_s) is calculated that is subsequently applied in the age equation, where , and P_P and P_S are the arithmetic mean of the weighted mean ²³⁸U/⁴³Ca ratios of each analysis of the zeta material (e.g., a Durango shard) from the primary and secondary LA-ICP-MS sessions, respectively. The zeta fractionation factor corrects for systematic (i.e., session-specific) fractionation in ²³⁸U/⁴³Ca ratios and can differ by up to 10% between sessions. This is related to variations in the analytical conditions (primarily related to the tuning of the ICP-MS instrument) between the primary LA-ICP-MS session and a given secondary LA-ICP-MS session. Applying the zeta fractionation, X_s, to the measured U/Ca ratio of each zeta analysis in a secondary LA-ICP-MS session means it can be directly compared with the corresponding values of the primary LA-ICP-MS session. This ratio should approximate unity for every zeta analysis, and if it shows systematic variation during a session, it is likely that the drift correction was inappropriate. If the drift correction was appropriate, then the session-specific zeta fractionation factor (X_s) is applied to all of the unknown grains analyzed during the secondary LA-ICP-MS session to generate fractionation-corrected and depth-weighted U/Ca ratios.

LA-ICP-MS–based fission-track dating is gradually overtaking the EDM in popularity due to its rapidity, the additional (trace-element and U-Pb age) information that can be acquired, elimination of the need for neutron irradiation, expediting turnaround times, mitigating radiation safety risks, avoiding administrative hassle, and negating the need to find one of the decreasing numbers of suitable reactors. In addition, mica etching requires HF, which represents a potential safety hazard.

One disadvantage of the LA-ICP-MS method is that it is a destructive technique on a microscopic scale, destroying part of each grain during analysis so that later inspection of the counting results is not possible. However, this shortcoming can be readily overcome if image sets are captured before ablation (i.e., the autonomous image-based approach described in Section 4.1), thus creating a permanent digital record of each grain (e.g., Gleadow et al., 2019).

The practical disadvantages of the EDM are offset by its appealing mathematical properties. The EDM is predicated on a so-called matched experimental design (Galbraith, 2010), which is based on the ratio of two Poisson variables, N_s and N_i. As a consequence, all of the statistical models can be fitted directly using the raw count data. Neutron dosimetry-based methods have the advantage that they are all based on a ratio of two-track density measurements, which have similar properties and are therefore likely to be subject to similar errors that ideally should cancel out.

In contrast, LA-ICP-MS–based fission-track data combine a Poisson-distributed numerator (N_s) with a denominator ([U] or [U/Ca]) that is assumed to be lognormal. The resulting ratios (or, equivalently, ages) are then approximated by a lognormal distribution in the different models outlined above. This approximation may not always be justified, especially when N_s is low.

In the experience of the authors, LA-ICP-MS zeta-based single-grain ages are generally more dispersed than EDM single-grain ages, both in terms of the statistical significance and absolute magnitude of the dispersion. The higher statistical significance of overdispersion reflects the higher precision of LA-ICP-MS data sets. This is caused by the low (percent-level) analytical uncertainty of mass spectrometric U-concentration measurements, as opposed to the higher counting uncertainties of induced fission-track densities using the EDM. The cause of the higher magnitude of LA-ICP-MS–based age dispersion is less clear, but likely reflects violations of the lognormal parametric assumption underlying the modified statistical models of Section 3.3. Additional dispersion is also likely to result from lateral zoning in U concentration between the spontaneous tracks counted in the area and the spot analyzed, although such errors are likely to be random, with minimal effect on the central age calculation.

The two analytical approaches are most likely to disagree in young and/or U-poor samples, especially in samples containing zero track grains. Whereas the Poissonian statistics behind the EDM can easily handle zero values, this is not the case for the LA-ICP-MS approach, whose lognormal assumptions do not allow zeros. Although Vermeesch (2017) shows that this problem can be circumvented by converting LA-ICP-MS data into “equivalent EDM” data sets, this can only be regarded as a stopgap solution whose accuracy remains unclear. Another problem for samples that are poor in spontaneous fission tracks is that it is impossible to determine, from the spontaneous fission tracks alone, whether the U concentration in a grain is uniformly distributed or strongly zoned. Thanks to its matched experimental design, such compositional zoning is more readily avoided using the EDM, whereas the LA-ICP-MS approach can only account for U zoning vertically (i.e., with depth) and not horizontally (i.e., along the surface of the grain mount). However, while induced track densities are also low in low-U samples, the accuracy of the mirror-image matching of areas counted for spontaneous and induced tracks is also compromised, which leads to additional errors if zoning is present.

Vermeesch (2017) proposed two strategies to mitigate the effects of U zoning. A first option is to only count fission tracks within the footprint of the laser spot that is used for the U measurement, thereby establishing the equivalent of a matched experimental design. The significant disadvantage of this approach is that it limits the track-counting area, resulting in a loss of precision. A second option, where possible, is to place multiple spots on some or all of the analyzed grains and use the dispersion of the U-concentration measurements as a measure of compositional heterogeneity. This approach does not hurt the precision but is less accurate than the matched experimental design option, because it assumes that the dispersion of the U measurements is the same for all grains in the sample. Any violation of this assumption will increase the dispersion of the data.

Ansberque et al. (2021) presented a LA-ICP-MS 2-D–mapping approach for apatites with very low spontaneous fission-track densities, with an analytical protocol that also facilitated the simultaneous acquisition of U-Pb age and trace-element data (Mn, Sr, La, Ce, Sm, Eu, Gd, and Lu). While the elemental mapping approach is ~10× slower than single-spot ablations, it is still faster than the EDM approach. Cogné and Gallagher (2021) undertook mathematical modeling employing random sampling of LA-ICP-MS fission-track data to compare the multispot and single-spot methods on apatite samples of known ages, variable zoning, and/or low-track densities. They showed that the multispot approach is robust for low-N_s and zoned samples, yielding both accurate and precise results without overdispersion. In contrast, single-spot measurements can induce overdispersion related to a relationship between single-grain age and U content. If single-spot analyses of apatites yield overdispersion when a single-age population may be expected, then, if the grain size is large enough, multispot analyses should be conducted to assess whether the excess dispersion is linked to undetected zoning and/or laser-spot misplacement, rather than to U-dependent annealing behavior. Available data suggest that zeta-derived EDM and LA-ICP-MS pooled ages from the same data set are generally in close agreement (Seiler et al., 2023; Fig. 6; see also Cogné et al., 2020), if they are based on similar assumptions regarding the etchable range of fission fragments. Table 2 provides a summary of the main pluses and minuses of EDM versus LA-ICP-MS.

Confined track lengths are measured routinely in apatite in combination with all of the different experimental methods for fission-track age dating. They are uniquely powerful recorders of thermal history, as each reflects track formation at a different point in time and annealing over the thermal conditions subsequently experienced, in aggregate that provides an integrated record of the time-temperature path of the host mineral grains. Confined tracks are classified based on the type of etchant pathway, predominantly as track-in-track (TINT; Lal et al., 1969) or track-in-cleavage (TINCLE; Bhandari et al., 1971; Figs. 1B, 1C, and 1E), with the former considered more reliable, as cleavages may provide pathways for geologic fluids to pre-etch or otherwise fix tracks and halt their annealing (Jonckheere and Wagner, 2000; Barbarand et al., 2003b).

When measuring confined tracks, only tracks that appear “fully etched” (Laslett et al., 1984) or “sufficiently etched” (Ketcham and Tamer, 2021) should be included. No formal definition exists for these, and if there were one, it would be complicated by the effects of anisotropy and variation in solubility on apatite. Generally speaking, a sufficiently etched track will have clear and easily distinguishable tips, and some indication that the along-track etch rate has slowed down to near-bulk etch rates, such as slow etching along internal crystallographic planes. Because this entails a judgement call and depends, in part, on training, it is highly recommended to calibrate track-length measurements for standards such as Durango and Fish Canyon against published studies. In particular, when using thermal history modeling, calibrating one's standards against the measured values that underlie the annealing models, and including a correction factor, is necessary to maximize the reliability of the derived histories (Ketcham et al., 2015, 2018).

The utility of track lengths in distinguishing among various thermal histories (e.g., fast cooling, slow cooling, and reheating) was first documented by Gleadow et al. (1986), and the advent of thermal history reconstruction software (Green et al., 1989; Issler, 1996a; Ketcham et al., 2000; Ketcham, 2005; Gallagher, 1995, 2012) made these interpretations more quantitative and detailed. In young or low-U apatite, ²⁵²Cf irradiation can be used to both increase the number of measurable confined lengths (Donelick and Miller, 1991) and reduce observer bias (Ketcham, 2005).

Due to the strong anisotropy of annealing, the resolution of track-length data can be increased by accounting for track angle using c-axis projection (Donelick et al., 1999; Ketcham, 2003; Ketcham et al., 2007a; see also Fig. 1B). C-axis projection removes anisotropic effects by normalizing track lengths for angle, correcting for analytical bias that may arise from differences in the orientation of measured confined tracks and thus providing more consistent thermal history reconstructions (Ketcham et al., 2009, 2018). The optimal c-axis projection model varies somewhat from analyst to analyst (Ketcham et al., 2007a, 2018), likely due to subtly different criteria used by analysts to measure tracks at different angles (Ketcham and Tamer, 2021).

With the advent of image-based data acquisition, there is renewed interest in attempting to utilize the length of tracks that intersect the grain surface, or semi-tracks. It is possible that their far greater numbers compared to confined tracks for a given sample could significantly improve statistical precision while still providing robust estimates of thermal histories (Laslett and Galbraith, 1996), though further theoretical development is required.

Fission tracks anneal as a function of time and temperature, with additional influential factors including crystallographic orientation (Green and Durrani, 1977; Donelick, 1991; Donelick et al., 1999; Ketcham, 2003); composition, e.g., in apatite, anions (Cl, F, and OH), and cation substitutions (e.g., REE, Mn, Sr, Fe, and Si; Green et al., 1985; Crowley et al., 1991; Carlson et al., 1999; Barbarand et al., 2003a); and, at least in zircon, radiation damage (e.g., Garver and Kamp, 2002; Rahn et al., 2004). Annealing can be generally viewed as a reconstruction of the crystal lattice to a sufficient degree that etching rates are no longer significantly enhanced, though the specifics of the lattice-scale mechanisms remain an area of active research. The earlier stages of annealing involve relatively gradual shortening from the track tips, whereas later stages are accelerated, with respect to temperature rise, and can include the formation of unetchable gaps (Green et al., 1986). Annealing rates are determined using laboratory experiments, where either spontaneous (also called fossil) tracks, or tracks created by reactor-induced fission of ²³⁵U, are heated for hours to months at high temperatures. Track lengths are measured and fitted to an annealing equation to extrapolate the resulting annealing trends to geological time scales; Ketcham (2019) provides a recent summary.

Although a quantitative physical understanding of fission-track annealing remains a desirable goal, the annealing equations that have proven most successful thus far in describing and extrapolating fission-track annealing are empirical (e.g., Laslett et al., 1987; Ketcham et al., 1999, 2007b). They are best visualized as a set of contours of constant annealing on an Arrhenius-type plot of the log of annealing time versus inverse temperature, as shown in Figure 7. The “fanning Arrhenius” model form (Laslett et al., 1987) posits a fanning set of iso-annealing lines and implies that the activation energy required for annealing increases with increasing extent of annealing. Subsequently, Ketcham et al. (1999) found that a “fanning curvilinear” form (first proposed by Crowley et al., 1991) did a better job of matching both high-temperature and low-temperature geological benchmarks for F-apatite (Fig. 7). The same appears to hold true for zircon fission tracks (Guedes et al., 2013; Ketcham, 2019). The curvature implies that the activation energy itself varies with temperature (Guedes et al., 2022) and can emerge in a variety of circumstances, such as when multiple individual reactions or processes contribute to a phase transition (Hulett, 1964; Carvalho-Silva et al., 2019). An associated simplification underlying current characterizations of fission-track annealing is the principle of equivalent time (Goswami et al., 1984; Duddy et al., 1988), which posits that the annealing rate of a fission track depends only on its length, and not its prior thermal history, which allows straightforward calculation of annealing along complex time-temperature paths. This assumption has come into question from a number of perspectives (Wauschkuhn et al., 2015; Ketcham and Tamer, 2021; Rufino et al., 2023), further emphasizing the need to develop a better physical understanding of the annealing process.

A fission-track age represents an integrated record of a mineral grain's thermal history starting from the time it last passed into a thermal regime in which fission tracks are retained. As such, a simplistic interpretation of an age as corresponding to an event, such as the time at which the sample cooled to below the system closure temperature, is only true in a limited number of cases (see Section 2). Instead, data are best interpreted in terms of the ways they elucidate the entire thermal history, and often with the aid of computational tools.

Along with the closure temperature (T_c; Dodson, 1973), a key concept for the interpretation of all fission-track data is the partial annealing zone (PAZ; Wagner, 1972; Wagner et al., 1989), a range of temperatures over which fission tracks anneal relatively rapidly. The canonical PAZ for fission tracks in apatite is ~60–120 °C (Gleadow, 1981), although this was a roughly calibrated interval that applies mainly to near-fluorapatite. The upper boundary of the PAZ corresponds to the total annealing temperature (Issler, 1996b), above which no fission tracks are retained at geological time scales. The lower boundary is defined by an inflection in down-well ages, or a strong deceleration, but not necessarily a halt in annealing rate, as fission tracks in apatite anneal even at Earth-surface temperatures (Donelick et al., 1990; Spiegel et al., 2007; Tamer and Ketcham, 2020; Vrolijk et al., 1992). The PAZ is best considered as the zone of maximum sensitivity of a fission-track system, but aspects of the thermal history after cooling below the PAZ can be extracted using confined lengths if they are well calibrated (e.g., Ketcham et al., 2018). The T_c is in the upper part of the PAZ. For near-end member F-apatite, T_c can range from 95 °C to 107 °C, according to experimental data for several different compositional varieties (Ketcham et al., 2007b).

On geologic time scales (10⁶–10⁷ yr), partial annealing and retention zones for the zircon fission-track system are approximated to be ~200–350 °C (Tagami et al., 1998; Rahn et al., 2004; Yamada et al., 2007; Ketcham, 2019). Increasing amounts of radiation damage are thought to lower annealing resistance in zircon (Garver and Kamp, 2002; Rahn et al., 2004). T_c estimates are 210–240 °C from field-based studies (Bernet, 2009), and somewhat higher (~280 °C) based on laboratory annealing experiments, with the difference possibly being that the latter utilized zircons with relatively low levels of radiation damage (Ketcham, 2019).

An important consideration often missed in routine usage of the PAZ and T_c is that annealing kinetics for apatite can vary greatly, as outlined in Section 4.5, and significantly impact annealing rates (see Fig. 8). Thus, it is risky to assign a strict thermal interpretation to AFT data without accounting for kinetics. The importance of Cl substitution in the anion site of apatite was documented in numerous field settings (e.g., Green and Duddy, 2012). The relative prevalence of other substitutions in sufficient magnitude to influence annealing kinetics is less well understood, as they are still seldom analyzed. Recent work in northwestern Canada (Powell et al., 2018) found that Cl alone was insufficient to explain age patterns, but a more complete accounting for kinetics that includes both OH content and cation substitutions (e.g., Carlson et al., 1999; Equation 6) does a better job. Fission-track etch pits can also be used to infer kinetics (Burtner et al., 1994) because solubility rises with Cl and OH content, although measuring D_par with sufficient reliability has proven unexpectedly difficult community-wide (Sobel and Seward, 2010; Ketcham et al., 2015, 2018).

Necessary preparatory steps for data interpretation include using the chi-squared test and radial plots to inspect whether the data contain multiple age and/or kinetic populations (see Fig. 3) and plotting single-grain ages and lengths against kinetic indicators, if measured. Detrital samples can easily feature multiple inherited components of different ages, and kinetic variation can occur in both sedimentary and crystalline rocks (e.g., O'Sullivan and Parrish, 1995). Kinetic variation will lead to a correlation between age and kinetic indicators, which will often match, but sometimes be clouded by, multiple inheritance; additional chemical data, such as REE profiles, or U-Pb age, or color or grain morphology, can be useful in disentangling the two factors (e.g., Chew and Donelick, 2012; Garver and Kamp, 2002).

A data set with multiple kinetic populations is advantageous because it is more information-rich, with multiple PAZs spanning a larger range of temperatures, with the different populations constraining each other because they need to be mutually consistent across their shared thermal history. Conversely, if a data set with multiple age populations is interpreted as a single population, the result runs the danger of being meaningless.

For sedimentary samples, or crystalline samples in close association with unconformities, a useful step is to compare the fission-track age to the stratigraphic age. If the fission-track age far predates the stratigraphic age, it can be considered non-reset postdeposition, or perhaps only lightly reset, and the data will contain information on the predepositional history, possibly including multiple inherited components. If the fission-track age is significantly younger than the depositional age, it was likely reset after deposition; only information on the postdepositional history (e.g., exhumation after burial) is present. If the fission-track and stratigraphic ages are similar (roughly within 20% of each other), it is likely that the fission tracks were partially reset and contain relatively precise temperature information. In samples with multiple kinetic populations, different populations may reflect different degrees of resetting, from unreset to fully reset.

Although the most detail can be extracted from confined length distributions using thermal history inversion methods, the basic patterns often can be discerned by eye. Broadly, a confined length distribution reflects the relative proportions of tracks that experienced temperatures in the PAZ versus those that did not, modified by observational biasing, as longer confined tracks are more likely to be intersected and measured than short tracks (Laslett et al., 1982). Some general patterns, modeled after the canonical classification set forth by Gleadow et al. (1986), are shown in Figure 9, for both nonprojected and c-axis–projected distributions. For fast cooling to below the PAZ, all track lengths are long. Slow cooling also features predominantly long track lengths, but with a short-track tail that reflects more time spent in the PAZ. Reheating paths can produce a bimodal pattern, reflecting tracks that formed before versus after peak burial. Extended residence time in the PAZ relative to subsequent cooling results in primarily short track lengths. In all cases, the modeled c-axis–projected tracks show the same pattern: while compressed over a shorter length interval, they have more distinct peaks due to the removal of dispersion caused by annealing anisotropy, resulting in overall better temperature resolution (Ketcham, 2003).

Often, fission-track data are interpreted via thermal history modeling software, predominantly using HeFTy (Ketcham, 2005) or QTQt (Gallagher, 2012). All such programs share the general approach of posing a series of candidate time-temperature paths and selecting those that best predict the observed age and length distribution, as well as any other data being modeled—e.g., (U-Th)/He ages and vitrinite reflectance Ro%—and any geological or other user-imposed time-temperature constraints. HeFTy and related programs, such as AFTINV (Issler, 1996a), use a frequentist approach, posing a null hypothesis that the model is correct, and evaluating the probability that the data could be observed given the model and measurement uncertainties. A Monte Carlo or iterative converging scheme is used to pose time-temperature paths, with the user specifying an allowed degree of complexity over different intervals of the history. QTQt uses Bayesian statistics, which allow the user to set up a series of priors, including annealing kinetic indicators. The Markov chain Monte Carlo method is used to search the solution space and includes a penalty for complex paths. Thus, it seeks and highlights the simplest set of solutions that reproduces the data. An important distinction is that if no paths fit the data to within uncertainties, HeFTy will return no result, whereas QTQt will return one regardless, requiring the user to evaluate whether the outcome explains the data acceptably by visually inspecting observed versus predicted age and confined track-distribution plots.

Both the frequentist and Bayesian approaches have their strengths and weaknesses (e.g., Vermeesch and Tian, 2014). However, similar advice applies for utilizing them and reporting their results. First, no program should be used as a black box; it is incumbent upon the user to understand what the program is doing, and all program inputs should be intentional and documented. Second, the geological context of the sample(s) being modeled should be understood, and only time-temperature constraints supported by robust independent evidence (e.g., stratigraphic age) should be enforced.

The reputation of thermal history modeling using apatite fission-track data was tarnished somewhat in the 1990s by artifacts in first-generation models, which utilized the annealing model constructed by Laslett et al. (1987) to describe experimental data by Green et al. (1986) for Durango apatite. The artifacts arose because the annealing model predicted insufficient annealing at near-surface conditions, requiring the inversion algorithm to maintain artificially high temperatures to achieve sufficient shortening of track lengths, which created an apparent late cooling event. The problem was remedied in subsequent annealing models by Ketcham et al. (1999, 2007b), which are based on annealing experiments on a wider range of apatites performed by Carlson et al. (1999) and Barbarand et al. (2003b) and use a different fanning curvilinear annealing equation. Durango apatite is more resistant to annealing than more common F-apatite, which generated part of the artifact, with the different annealing equation accounting for the remainder. With these corrections, the artifact was eliminated, and the <60 °C portions of thermal histories (e.g., Shorten and Fitzgerald, 2021) are much more reliable, although attention to length-measurement calibration is important for maximizing confidence (Ketcham et al., 2015, 2018).

Fission-track data can also be interpreted thermo-kinematically, using software such as Pecube (Braun et al., 2012; McQuarrie and Ehlers, 2015) or FetKin (Almendral et al., 2015; Mora et al., 2015) that creates a 2-D or 3-D thermal model that simulates a situation of interest (e.g., tectonic exhumation or topographic development) and uses the time-temperature paths calculated for various locations in the model space to predict ages and mean lengths. In some implementations, models may iterate over various parameters (e.g., fault timing and slip rate) to attempt to better fit the data (e.g., Fox et al., 2014; Valla et al., 2010). These methods tighten the link between the data and the geological interpretation, but at some cost in detail, particularly in the case of fission-track data, where current implementations fit mean track lengths but not the complete length distributions. They are also limited in that any inherited history (i.e., fission tracks formed before the start of the thermal model) can only be very approximately included. However, if the problem is well-posed to avoid being affected by such limitations, these methods enable geological questions to be addressed much more directly.

Fission-track analysis of detrital apatite and zircon grains from modern river, beach, or deltaic sediments and/or ancient sandstone is a widely used technique for studying (1) the exhumation history of orogenic mountain belts, such as the Himalayas (e.g., Cerveny et al., 1988), European Alps (e.g., Bernet et al., 2009), and the Andes (e.g., Bermúdez et al., 2013; Coutand et al., 2006); (2) provenance analysis (e.g., Hurford and Carter, 1991), particularly in combination with single-grain, fission-track/U-Pb double dating (e.g., Carter and Moss, 1999) and trace-element geochemical analysis of double-dated grains (Hülscher et al., 2021); and (3) the thermal evolution of sedimentary basins (e.g., Naeser et al., 1989). Samples of detrital apatite and zircon from modern river and beach sediments were demonstrated to provide a relatively good overview of the spectrum of fission-track cooling ages of bedrock exposed in the drainage area. Therefore, a fairly reliable first-order overview on exhumation rates and sediment provenance can be obtained (e.g., Bermúdez et al., 2013). Nonetheless, potential lithological, exhumation rate, and transport biases must be considered (Malusà et al., 2013; Naylor et al., 2015). For ancient clastic sediments, when interpreting fission-track data from detrital sedimentary rocks, we must also consider the influence of postdepositional burial heating, and when working with detrital samples from outcrops, the cooling history during basin inversion must be considered.

Because the detrital apatite and zircon grains in a sediment sample may be derived from a wide range of source rocks with more or less complex individual thermal histories, it is best to date at least 100 or more individual grains per sample to cover the whole grain-age spectrum (Vermeesch, 2004). In many studies that is not possible, so if the objective is to determine sediment provenance and the exhumation history of the source area, as many grains as possible should be dated. To determine the thermal history of sedimentary basins and to study partial or full postdepositional annealing of fission-tracks, at least 20 grains or more should be dated, and with up to 100 track-length measurements, if possible, obtained on many more grains from the same sample.

For detrital samples, use of the central age is recommended (Galbraith and Laslett, 1993), as grain-age distributions of detrital samples tend to be overdispersed (see Section 3.3; Fig. 3). This may not be the case, however, when detrital grains were completely thermally reset with respect to the temperature sensitivity of the applied fission-track system postdeposition, which is reflected in high P(χ²) and low age dispersion (<20%) values. The central age of a detrital sample can be used to obtain a crude estimate of drainage basin mean erosion rates (e.g., Bernet et al., 2009), or compared with the depositional age to evaluate the degree of partial annealing due to postdepositional burial heating (see below). Detrital apatite or zircon fission-track data presented in research articles need to include all commonly reported data in fission-track studies, and this will vary depending on whether the external detector method or the LA-ICP-MS method is used (see Section 6).

Because of the commonly wide range of observed grain ages in detrital samples, it is possible to decompose the observed grain-age distribution into major population ages by peak fitting (Galbraith and Green, 1990), as discussed below. The use of best-fit peak ages is not obligatory in detrital thermochronological studies, given that in many cases, as stated previously, the meaning of particular peak ages is uncertain, as a single peak age normally does not simply relate to a specific cooling or tectonic event. Such interpretations need to be based on a series of analyses from the stratigraphic record to determine static or moving age peaks that can be tied into the geological evolution of a study area (Bernet and Garver, 2005). If peak-fitting is used, then the peaks should be grouped consecutively based on peak age, with indication of the percentage of grains belonging to each peak, to compare the results from different samples more easily. A peak should contain at least 5% of the dated grains for it to have some significance. For peak-fitting of detrital fission-track data, it is important to inspect the results, as they should lead to geologically meaningful interpretations. This inspection includes checking how much individual peaks overlap or if they are well separated. Also, the number of grains analyzed is important, as a larger number of dated grains permits more peaks to be fitted (Vermeesch, 2019), even though they may not convey any geological information.

Fission-track data obtained from detrital sediments, whether they are modern river or beach sediments or ancient sandstone, are notoriously complicated, because they may represent a continuous or discrete mixture of cooling ages. Within the same drainage area, it is possible that bedrock with preorogenic and postorogenic cooling ages as well as partially reset cooling ages is present, with both providing detrital grains to analyzed sediments and sedimentary rocks. The challenge is to interpret the observed grain-age distribution of a detrital sample in a geologically meaningful way. The mixing of non-reset and partially reset grains from source rocks with grains carrying orogenic (e.g., related to postmetamorphic exhumation) apparent cooling ages commonly results in an overdispersed, continuous mixture of grain ages in both modern and ancient sediments (Figs. 3 and 4). Furthermore, because of burial heating, detrital apatite may be particularly affected after deposition due to its relatively low-temperature PAZ.

A common approach for analyzing detrital grain-age distributions is to identify different age groups or clusters that may be indicative of sediment provenance, peak exhumation rates, or the timing of tectonic events. For this purpose, different peak-fitting techniques were proposed. The most widely used software packages for this type of analysis over the past 10–20 years were the BINOMFIT program of Brandon (see Ehlers et al., 2005) and the RadialPlotter program of Vermeesch (2009). Both programs provide an automatic fit option for determining best peak-fit solutions, but also allow manual selection of the number of peaks that should be fitted to the data set. The second option may be useful if automatic fitting provides a solution with considerable peak overlap. Geological constraints should have the highest priority for constraining models and their interpretation. Depending on the objective of the study (e.g., provenance analysis and comparison to bedrock data sets), peak fitting may be useful, but the interpretation of individual peaks needs to be done with caution, as not all peaks necessarily reflect a tectonic or rapid cooling event. For tracing rapid cooling or exhumation, it is safer to concentrate just on the youngest age peak or the so-called minimum age (Galbraith and Laslett, 1993; Section 3.3; Fig. 4). The minimum age is an estimate of the first coherent age component in an observed grain-age distribution. This approach avoids the sample size bias to younger ages by peak-fitting (Vermeesch, 2019).

In this respect, it is also necessary to check whether the grains that belong to the minimum age group in a detrital sample were affected by postdepositional partial annealing; comparing the difference between the stratigraphic age (based on independent stratigraphic information) and the minimum age provides a first indication. In the case where the minimum age peak predates but does not overlap with the stratigraphic age within statistical uncertainties, no (significant) postdepositional partial annealing can be expected. If the minimum age or the first peak overlaps with the stratigraphic age, one needs to check whether rapidly cooled volcanic apatite or zircon grains may be present. This could be confirmed or rejected by double dating with U-Pb rim ages of the fission-track–dated grains (e.g., Jourdan et al., 2013; Mark et al., 2016). If no volcanic grains are present, then one needs to take the stratigraphic position, possible maximum burial depth, and basin thermal gradients into account to distinguish between rapid source-area exhumation or partial annealing due to deposition and burial. For detrital apatite fission-track analysis, the additional measurement of horizontal confined tracks is a commonly used technique (see Section 4.5), but where possible these should be acquired from the same grains from which age determinations were made.

For detrital apatites and zircons from ancient sedimentary rocks that were unaffected by partial annealing during burial, the lag-time concept can be applied to estimate long-term sediment source-area exhumation rates. In this context, the lag-time is defined as the difference between the fission-track peak or minimum age and the age of sediment deposition in the basin (e.g., Garver et al., 1999; Bernet and Garver, 2005). To apply the lag-time concept, the stratigraphic age needs to be independently constrained through magnetostratigraphic and/or biostratigraphic data.

Finally, if the minimum age or the first peak are clearly younger than the stratigraphic age, then significant postdepositional partial annealing is commonly the best interpretation. Certain basin studies showed that samples with minimum ages, or first-peak ages, or even central ages collected from an inverted basin may line up around a common age with increasing paleo-depth. Such samples may indicate the timing of basin inversion if the rocks cooled from T_c to the surface in a single event, which means that exhumation was significant and cooling rapid (e.g., Sobel and Dumitru, 1997; van der Beek et al., 2006). Older age peaks determined in partially annealed basin samples may still have relatively long lag times but reflect a mixed cooling signal from within the PAZ and the inherited cooling signal from source-area exhumation. Such older peak ages must be interpreted with care and are not recommended for calculating exhumation rates.

Fission-track analysis of minerals from bedrock samples was also widely applied to geologically challenging areas, such as cratons or passive continental margins, which commonly offer limited or no syn- or post-tectonic geological record with which to inform thermochronological interpretations (e.g., Gallagher and Brown, 1997; Gallagher et al., 1998; Persano et al., 2005; Wildman et al., 2015; Enkelmann and Garver, 2016; Kohn and Gleadow, 2019). For example, as indicated previously, the use of apatite fission-track analysis to quantitatively reconstruct thermal histories via the inversion of track-length data proved to be very effective in such applications, as it allows one to decipher and constrain the timing, rate, and magnitude of exhumation, and in some cases, burial events, even in the absence of preserved geological strata that directly record that event. However, one must acknowledge that the fission-track method has an intrinsic limitation in that the modeled thermal histories reproducing measured fission-track data can be nonunique. Therefore, results from the inversion of fission-track data largely depend on predefined time-temperature constraints, which in the case of cratons or passive continental margins are often only inferred given the frequent absence of independent time-temperature constraints that can be incorporated into models.

The Archean Pilbara Craton, northwestern Australia, provides an example where two different thermal histories were derived from modeling just an AFT data set (Fig. 10A) and then another where additional constraints were applied to remodel the initial data (Fig. 10B). The addition of a rare geological constraint in this setting, together with coexisting apatite (U-Th-Sm)/He data, reveal a pronounced Paleozoic cooling episode, which is also recorded from other localities across the Pilbara Craton, suggesting the removal of a kilometer-scale sedimentary blanket that is thought to have been triggered by the widespread intraplate Paleozoic Alice Springs Orogeny (Morón et al., 2020). This example underscores the importance of using well-defined time-temperature constraints for thermal history modeling, where possible.

The passive continental margin of West Greenland presents an example where there are contrasting views on the tectono-thermal evolution of the region due to differences in thermochronological data interpretations (e.g., Japsen et al., 2006; Redfield, 2010; Jess et al., 2018; Danišík and Kirkland, 2023, and references therein). West Greenland was subject to Cretaceous rifting and Paleocene breakup from Canada (e.g., Oakey and Chalmers, 2012), and several works proposed that the present-day topography of central West Greenland reflects several episodes of post-breakup uplift and erosion. This so-called polycyclic model of passive continental margin formation is built on several independent lines of evidence, which include thermal history modeling of apatite fission-track data. For example, for a basement sample (code: GC891-27) with an AFT age of 113 ± 9 Ma and mean track length (MTL) of 11.5 ± 2.3 μm, Green et al. (2003) modeled an episodic thermal evolution with two heating episodes and two cooling episodes with onset at 160–115 Ma and at 40–0 Ma (Fig. 11). This result was achieved by applying time-temperature constraints inferred from the geological record in adjacent or remote areas. In contrast, Jess et al. (2018) remodeled the same data set without applying any predefined constraints, arguing for an absence of the geological constraints directly in the investigated location, and suggested a history of monotonous cooling from the apatite PAZ (60–120 °C; Wagner and Van den haute, 1992) since ca. 180 Ma until the present (Fig. 11). These modeling results imply that the elevated topography of central West Greenland formed in the Mesozoic and remained elevated due to the isostatic response to protracted erosion (Redfield, 2010; Jess et al., 2018). This example from West Greenland demonstrates that the strategy for inversion of fission-track data is critical for interpretational outcome. There are numerous reasonable approaches to modeling the fission-track data, and their applications depend on the type of data, interpretation of geological context, availability of additional information from other geo-thermochronological methods, or research hypothesis being tested by the thermal modeling.

The challenge in making recommendations for (meta-)data reporting lies in designing schemas that are standardized while remaining flexible enough to accommodate fission-track data produced via a range of rapidly advancing analytical techniques. Depending on the specific aims, and the procedures and analytical instruments employed in a given study, the types of (meta-)data reported will naturally vary. Therefore, such schema must be pliable enough to fit both analog and digital (manual and image-based) counting methods and confined track measurements, a range of age-determination methodologies and associated age equations (see Sections 2 and 3), and a variety of analytical techniques that are increasingly used to determine in situ parent ²³⁸U concentrations and geochemistry (e.g., LA-ICP-MS and electron probe microanalysis [EPMA]).

Nevertheless, certain common data reporting principals exist that can guide FAIR fission-track data reporting. Regardless of the specific analytical techniques employed, the responsible reporting of fission-track (meta-)data and associated analyses should:

• Include any geosample information needed to place the data in a geological context.

• Allow readers to independently repeat experimentation, including recalculating ages and performing thermal history modeling.

• Enable the reader to independently assess data quality through consistent reporting of uncertainties, calibration parameters, totals, and secondary reference material results.

• Make detailed (meta-)data available in a tabulated format to enable easy extraction and reuse by both humans and machines, in addition to reporting data in the main text.

Table 3 provides a checklist to guide authors and reviewers in how to best present their fission-track and associated results following FAIR data principals. In the subsequent sections, recommendations for minimum required and highly recommended (meta-)data reporting for geosamples, fission-track counting and length measurements, mineral composition analyses and kinetic parameters, fission-track age determinations, and thermal history models are presented. While the number of (meta-)data we recommend be required for FAIR fission-track data reporting is large, the increase in the use of field geology applications, digital fission-track analytical software (TrackWorks and FastTracks; Gleadow et al., 2009), and data reduction software (e.g., Iolite; Paton et al., 2011) makes the capturing and reporting of these data increasingly automated and streamlined.

Increasingly, scientific journals require authors to present detailed (meta-)data in open access data repositories. Data platforms that allow users to share, access, explore, and interrogate geochemical data from around the globe in a 4-D (3-D + time) geospatial context, such as AusGeochem (Boone et al., 2022), provide a powerful tool for FAIR thermochronological data dissemination. Its bespoke relational data schemas for fission-track and (U-Th)/He thermochronological analyses (Boone et al., 2023) match the data reporting guidelines presented below and enable a range of real-time data plotting and interpolation tools.

Fission-track data are only as good as the geospatial and geological contexts in which they can be placed. If these are not robust, the thermal histories that they record are often cryptic and their geoscientific implications difficult to interpret. Therefore, it is critical for detailed geosample metadata to be reported along with their corresponding fission-track data.

The most critical geosample metadata (Table 4) include the sample identification (ID), the latitude-longitude of the sample along with the geodetic datum used, elevation and vertical datum, sample lithology, mineral type, sample kind and method, and sample depth in the case of core, well, or dredge samples. We also strongly encourage data producers to mint International Geo Sample Numbers (IGSNs; Lehnert et al., 2011; Klump et al., 2021) for their samples in addition to the sample ID generally given to a sample in the field. The use of globally unique and persistent IGSN identifiers for physical samples prevents unintentional sample ID duplication in global repositories and enables samples to be discoverable on the internet (Klump et al., 2021).

Additional recommended geosample metadata include latitude/longitude precision; detailed information about the sampling location, stratigraphic unit name, and chronostratigraphic unit age; any independent constraints for crystallographic or metamorphic age of the rock and/or mineral sample; sample collector, collection date, and last known archive of the sample; and any prior publications associated with the sample.

Detailed fission-track counting data should be reported on a per-grain basis, enabling single-grain, fission-track ages to be independently recalculated and modeled. Minimum single-grain counting data to report include where and by whom the analyses were performed, the detailed analytical methodology, the etching conditions, all relevant counting statistics, average per-grain D_par values, and, in the case of EDM results, the dosimeter used (Table 5). Grain IDs should also be reported so that fission-track counting results can be related to any determinations of semi- and confined fission-track length (Section 6.3); corresponding mineral compositional and kinetic parameter data (Section 6.4); and calculated single-grain, fission-track ages (Section 6.5) from the same grains. Analytical metadata pertaining to the analyst, laboratory, etching conditions, microscope, objective, and the fission-track counting method (manual or image-based) should also be presented. In the case of digital fission-track count data, additional required metadata include the analytical software and version, and the fission-track analytical algorithm employed. Reporting sample mount IDs, the date of analysis, the counting area in square centimeters, and D_per values (if possible) is also recommended for completeness.

In addition to presenting mean track lengths, standard errors, standard deviations, and the number of tracks measured per sample, it is critically important that confined track measurements are reported on a per-track basis (Table 6), along with semi-tracks, if measured. Such detailed per-track data are required to enable readers to independently perform thermal history modeling using published results. These should be accompanied by analytical metadata describing the methodology employed, such as the microscope and objective used, the etching conditions, and whether samples were irradiated with ²⁵²Cf or heavy ions to increase the number of confined track measurements (e.g., Donelick and Miller, 1991; Jonckheere et al., 2007). If track measurements were performed digitally, the software and version used should be disclosed. In the case of digitally obtained fission-track measurements, other recommended data to report include track dip, azimuth, and corrected z-depth.

Increasingly, in situ mineral compositional data are obtained in tandem with digital fission-track analysis, both for determining the parent ²³⁸U content of grains analyzed for fission-track age calculations and for measuring granular and subgranular major, minor, and trace-element concentrations. Such compositional data are generally produced by some combination of LA-ICP-MS and EPMA. In addition to reporting these geochemical results, authors should also take care to report the associated LA-ICP-MS (Table 7) and EPMA (Table 8) analytical metadata required for readers and data users to be able to reuse and independently assess data quality.

Any compositional and kinetic indicators should also be presented, along with the equations and parameter values used for calculation (Table 9). These may include r_mr0 and the associated κ parameter (Carlson et al., 1999; Ketcham et al., 1999), and transformations of r_mr0 determinations to make them more linear and intuitive (Issler et al., 2022), such as effective chlorine (eCl) or effective D_par (eD_par). Authors may also wish to report other compositional parameters, such as effective U content, a proxy for accumulated radiation damage that is sometimes used to assess the potential effect of α-radiation in enhancing fission-track annealing (e.g., McDannell et al., 2019).

Like confined track and semi-track data, detailed fission-track age determinations must be reported on a per-grain basis to ensure that numerical thermal history modeling can be independently performed (Table 10). These should be accompanied by statistical measurements of single-grain age variance, reported as dispersion, a chi-square test—P(χ²)—or preferably both.

Given the different fission-track age equations and approaches being used by fission-track laboratories around the globe (discussed in Sections 3 and 4), the equations, constants, parameters, and calibration factors used for age calculation must also be reported. Doing so will enable future work to equate regional fission-track data by recalculating ages in bulk, for example, using updated constants or similar R_s or q factors, a prospect that is becoming much more realistic with the advent of structured relational thermochronological databases (Boone et al., 2023).

In the case of detrital fission-track analysis, it is also recommended that authors report the ages, uncertainties, and the grain percentage for each best-fit peak, if determined. The program used to obtain peaks must be listed, and the peak-fitting approach used (automated or user-defined number of age peaks) should be described in the methodology.

Thermal history modeling of fission-track and other thermochronological data is now routine in studies that employ these methods. Nevertheless, large variations exist in the manner and degree of reporting of thermal history modeling (meta-)data. This prompted Flowers et al. (2015) to recommend a reporting protocol for thermochronological modeling. In that work, and the subsequent comment (Gallagher, 2016) and response to comment (Flowers et al., 2016), a framework for reporting modeling metadata and results was established to enable reviewers and readers to independently judge numerical simulation quality and validity. This framework is summarized in Table 11.

We emphasize that there is a strong need in scientific papers for authors to clearly explain the strategy for interpretation of fission-track data, and to report all data acquired, allowing readers to reproduce the results or develop alternative models (e.g., full single-grain age data and full track-length data, including kinetic indicator measurements such as D_par, Cl, or r_mr0 for apatite; see Tables 5–10 for more details). In addition, we recommend that publications should state the following: the objective of the modeling, e.g., why the modeling was undertaken, or what research hypotheses are being tested, and the rationale for the model setup. Detail to be provided in the rationale for the model setup includes the type of time-temperature constraints that were applied, particularly those from independent geological information, and why; setting parameters, such as the annealing model employed, initial track length, kinetic parameter used, and the statistical fitting criteria; and how well the resulting simulations reproduce the observed data. Many studies model other thermochronological data together with fission-track data, and in such cases the other data sets should also be appropriately documented (e.g., Flowers et al., 2023b).

With the number of fission-track analyses acquired globally rapidly growing, and the scale and scope of fission-track investigations expanding, the need for consistent and detailed fission-track data reporting is greater than ever. To further advance the burgeoning field of AI-powered, automated fission-track analysis (Nachtergaele and De Grave, 2021; Li et al., 2022)—for instance, ever larger data sets of 2-D and 3-D photomicrographs of fission tracks in mineral grains and micas—will be required, along with the machine-readable (meta-)data that describe the micron-scale geometries of the thousands to millions of individual fission tracks that they contain. At the same time, increasingly large thermochronological data sets are being interrogated on regional (e.g., Kohn et al., 2002, 2005; Boone et al., 2021) to global scales (Herman et al., 2013), demanding that comprehensive fission-track data be reported on a per-grain and even per-track basis to enable regional age recalculation and thermal history modeling and remodeling to be performed using the same parameters.

Fortunately, an ever-growing range of digital fission-track analysis (TrackWorks and FastTracks, Gleadow et al., 2009; AI-Track-tive, Nachtergaele and De Grave, 2021), LA-ICP-MS data reduction (e.g., Iolite, Paton et al., 2011), numerical modeling (e.g., HeFTy, Ketcham, 2005; QTQt, Gallagher, 2012), and geospatial data platform (AusGeochem, Boone et al., 2022, 2023) programs are available to help expedite the collation and dissemination of enormous volumes of detailed fission-track (meta-)data. Nevertheless, additional software development is required to streamline the data continuum among these disparate programs through mechanisms such as application programming interfaces (APIs), which have the potential to automatically capture data from the laboratory and store them in FAIR data repositories. In the interim, an unavoidable degree of manual labor is required of the thermochronological community to ensure that their detailed fission-track (and other thermochronological) data are Findable, Accessible, Interoperable, and Reusable.

Increasing levels of automation in fission-track analysis also include the promise of being able to capture additional track-measurement parameters, such as detailed 3-D anisotropic track dimensions and orientations (e.g., Jonckheere et al., 2019; Aslanian et al., 2021). New analytical technologies should also allow simultaneous multielement compositional analysis alongside the track measurements. Combined with a new generation of anisotropic etching, track annealing, and thermal history models, such additional parameters have the potential to further enhance the interpretations that can usefully be extracted, but at the expense of a further increase in the complexity of data sets. Indeed, such additional measurements could only be contemplated in any practical sense because of the emergence of automated data acquisition techniques. We anticipate, therefore, that the already obvious demand for more comprehensive and standardized approaches to data reporting will continue to increase.

The burden to enact and enforce this much-needed change in data reporting practice will ultimately fall to the thermochronological community itself, as both data providers and scientific reviewers. The benefits of doing so, however, are extremely exciting and have the potential to take the field of (fission-track) thermochronology to new horizons. The methodical archiving of detailed thermochronological (meta-)data in structured schemas will facilitate the rapid interrogation of intractably large data sets comprising hundreds of analyses produced by numerous laboratories from around the globe. This, in turn, can enable future integration of regional thermochronological data with other numerical geoscience techniques driving new scientific insights, such as bulk thermal history modeling, numerical landscape evolution modeling, and the placement of fission-track analyses in palinspastic and/or paleoclimatic reconstructions. Fission-track data reporting guidelines should remain an ongoing topic of community discussion as further developments occur.

The University of Melbourne Thermochronology Laboratory receives support from the AuScope program (www.auscope.org.au) of the Australian government's National Collaborative Research Infrastructure Strategy (NCRIS). The authors thank Willy Guenthner and an anonymous reviewer for their constructive comments and suggestions and Brad Singer for editorial handling.