gdtchron.aft

Module for forward modeling of the apatite fission track system.

This code follows the workflow from Ketcham (2005) and supports producing ages and distributions of c-axis projected lengths. This code uses the fanning curvilinear (FC) model from Ketcham et al. (1999).

Functions

`calc_aft_age`(r_final, tsteps[, rho_st])	Calculate AFT age based on present-day reduced lengths.
`calc_annealing`(r_initial, temperature, ...)	Calculate the annealing of fission tracks across a timestep.
`calc_l_dist`(r, dpar, tsteps[, constants, ...])	Calculate length distribution from reduced lengths.
`calc_weights`(r, tsteps[, lamb])	Calculate weights to apply to length distributions when summing together.
`combine_dists`(means, stdevs, w[, ...])	Calculate weighted sum of normal distributions.
`dpar_conversion`(r_mr, dpar[, constants])	Convert reduced lengths for one apatite to another apatite.
`f`(temperature, time_annealed[, constants])	Calculate f following Equation 4 from Ketcham et al. (1999).
`forward_model_aft`(temps, tsteps, dpar[, ...])	Conduct forward modeling of the apatite fission track system.
`g`(r[, constants])	Implement the length transform from Ketcham et al. (1999) (Equation 6).
`get_equiv_time`(r_initial, temperature[, ...])	Calculate time it takes to reach a reduced length at given temperature.
`get_l_stdev`(l_c)	Get standard deviation of length distribution based on its mean length.
`get_next_r`(temperature, time_annealed[, ...])	Calculate reduced lengths after annealing over a given time period.
`l_conversion`(r, dpar[, constants])	Convert reduced length r to mean c-axis projected length l_c.
`r_to_rho`(r)	Convert final reduced lengths to fission track densities.

gdtchron.aft.calc_aft_age(r_final, tsteps, rho_st=0.893)[source]

Calculate AFT age based on present-day reduced lengths.

This function uses Equations 13-14 from Ketcham (2005). Note that, if rho_st < 1, then an un-annealed apatite will have an age older than the earliest timestep. This behavior is reproduced in HeFTy when using both the Ketcham et al. (1999) and Ketcham et al. (2007) parameters.

Parameters:

r_final (NumPy array of floats) – Mean reduced lengths (unitless) of fission tracks produced at each point on the apatite’s time-temperature path. These values should already be converted to account for dpar variations.
tsteps (NumPy array of floats) – Array of times (yrs BP) that each set of fission tracks was produced at. This array should be in descending (i.e., chronological) order. The final time in this array should be 0. and should not have a corresponding reduced length in r_final.
rho_st (float) – Fission track density reduction in the age standard (default: 0.893, the value for the Durango apatite)

Returns:

AFT age (Myrs BP)

Return type:

float

gdtchron.aft.calc_annealing(r_initial, temperature, start, end, next_nan_index, constants={'alpha': -0.12327, 'beta': -11.988, 'c0': -19.844, 'c1': 0.38951, 'c2': -51.253, 'c3': -7.6423, 'l_intercept': 15.72, 'l_slope': 0.35, 'r_kappa_sum': 1.0})[source]

Calculate the annealing of fission tracks across a timestep.

Parameters:

r_initial (NumPy array of floats) – Initial mean reduced lengths (unitless) of fission tracks at start of timestep. The value at index 0 corresponds to fission tracks produced at the first timestep, the value at index 1 corresponds to fission tracks produced at the second timstep, and so on. np.nan should be stored at indices corresponding to fission tracks produced at the current timestep or at future timesteps.
temperature (float) – Temperature (K)
start (float) – Start time of timestep (yrs BP)
end (float) – End time of timestep (yrs BP)
next_nan_index (int) – First index of r_initial to have a value of np.nan. This index corresponds to fission tracks that will be produced at the current timestep.
constants (dict) – Dictionary of constants associated with annealing model being used (default: KETCHAM_99_FC)

Returns:

Updated mean reduced length(s) of fission tracks at the end of this timestep.

Return type:

NumPy array of floats

gdtchron.aft.calc_l_dist(r, dpar, tsteps, constants={'alpha': -0.12327, 'beta': -11.988, 'c0': -19.844, 'c1': 0.38951, 'c2': -51.253, 'c3': -7.6423, 'l_intercept': 15.72, 'l_slope': 0.35, 'r_kappa_sum': 1.0}, make_graph=False, l_num=100)[source]

Calculate length distribution from reduced lengths.

This function takes the reduced lengths of fission tracks produced at each timestep, calculates the mean and standard deviation of the length distributions for each timestep, and combines these distributions to find the overall fission track length distribution

Parameters:

r (NumPy array of floats with length n) – Reduced lengths of fission tracks produced at each timestep (unitless)
dpar (float) – Etch figure length (micrometers)
tsteps (NumPy array of floats with length n + 1) – Array of times (yrs BP) that each set of fission tracks was produced at. This array should be in descending (i.e., chronological) order. The first time is the start of the first timestep, last time is end of last timestep.
constants (dict) – Dictionary of constants associated with annealing model being used (default: KETCHAM_99_FC)
make_graph (bool) – Boolean indicating whether or not to return x and frequency series for plotting (default: False)
l_num (int) – Nonnegative integer specifying how many lengths to use for the frequency series and internal calculations (default: 100). A value of at least 50 is recommended for calculating the standard deviation and mean of the mixed distribution; larger values produce smoother curves when plotting the frequency series.

Returns:

mean (float) – Mean of mixed distribution
stdev (float) – Standard deviation of mixed distribution
l_c (optional) (NumPy array of floats) – Array of length values. Only returned if make_graph = True
freqs (optional) (NumPy array of floats) – Array of frequencies at which each x value is observed. Only returned if make_graph = True

gdtchron.aft.calc_weights(r, tsteps, lamb=0.0001551)[source]

Calculate weights to apply to length distributions when summing together.

This function uses Equations 12 and 13 in Ketcham (2005) to determine the weights to apply to the length distributions associated with each timestep. These weights are used when summing the length distributions together to create a mixed distribution.

Parameters:

r (NumPy array of floats with length n) – Mean reduced lengths (unitless) of fission tracks produced at each point on the apatite’s time-temperature path. These values should already be converted to account for dpar variations.
tsteps (NumPy array of floats with length n + 1, where n > 1) – Array of times (yrs BP) that each set of fission tracks was produced at. This array should be in descending (i.e., chronological) order. The first time is the start of the first timestep, last time is end of last timestep.
lamb (float) – Total U238 decay constant (lambda) (Myr^-1) (default: 1.551e-4)

Returns:

Weights for each timestep

Return type:

NumPy array of floats with length n

gdtchron.aft.combine_dists(means, stdevs, w, make_graph=False, x_num=100)[source]

Calculate weighted sum of normal distributions.

This funcion calculates a weighted sum of normal distributions based on their means and distributions. Disregards means of zero and negative values. Used to create mixed distribution of lengths

Parameters:

means (NumPy array of floats) – Means of each distribution
stdevs (NumPy array of floats) – Standard deviations of each distribution
w (NumPy array of floats) – Weights to apply to each distribution when combining
make_graph (bool) – Boolean indicating whether or not to return x and frequency series for plotting (default: False)
x_num (int) – Nonnegative integer specifying how many x values to use for the frequency series and internal calculations (default: 100). A value of at least 50 is recommended for calculating the standard deviation and mean of the mixed distribution; larger values produce smoother curves when plotting the frequency series.

Returns:

mean (float) – Mean of mixed distribution
stdev (float) – Standard deviation of mixed distribution
x (optional) (NumPy array of floats) – Array of x values. Only returned if make_graph = True
freqs (optional) (NumPy array of floats) – Array of frequencies at which each x value is observed. Only returned if make_graph = True

gdtchron.aft.dpar_conversion(r_mr, dpar, constants={'alpha': -0.12327, 'beta': -11.988, 'c0': -19.844, 'c1': 0.38951, 'c2': -51.253, 'c3': -7.6423, 'l_intercept': 15.72, 'l_slope': 0.35, 'r_kappa_sum': 1.0})[source]

Convert reduced lengths for one apatite to another apatite.

This function converts from the reduced lengths of a more resistant apatite to reduced lengths of a less resistant apatite following Equations 7, 8, and 9a of Ketcham (2005).

Parameters:

r_mr (NumPy array of floats) – Mean reduced lengths (unitless) of fission tracks for the more resitant apatite
dpar (float) – Etch figure length (micrometers) for the apatite
constants (dict) – Dictionary of constants associated with annealing model being used (default: KETCHAM_99_FC)

Returns:

Mean reduced lengths (unitless) of fission tracks for the more resitant apatite.

Return type:

NumPy array of floats

gdtchron.aft.f(temperature, time_annealed, constants={'alpha': -0.12327, 'beta': -11.988, 'c0': -19.844, 'c1': 0.38951, 'c2': -51.253, 'c3': -7.6423, 'l_intercept': 15.72, 'l_slope': 0.35, 'r_kappa_sum': 1.0})[source]

Calculate f following Equation 4 from Ketcham et al. (1999).

Parameters:

temperature (float) – Temperature (K)
time_annealed (float or NumPy array of floats) – How long the crystal annealed at a given temperature (s)
constants (dict) – Dictionary of constants associated with annealing model being used (default: KETCHAM_99_FC)

Returns:

Value(s) of f for each value of t

Return type:

float or NumPy array of floats

gdtchron.aft.forward_model_aft(temps, tsteps, dpar, constants={'alpha': -0.12327, 'beta': -11.988, 'c0': -19.844, 'c1': 0.38951, 'c2': -51.253, 'c3': -7.6423, 'l_intercept': 15.72, 'l_slope': 0.35, 'r_kappa_sum': 1.0}, get_lengths=False, make_graph=False, l_num=100)[source]

Conduct forward modeling of the apatite fission track system.

This function runs the forward model of the AFT system from Ketcham (2005) based on given time-temperature series. Calculates thermochronological age and (optionally) mean fission track length

Parameters:

temps (NumPy array of floats with length n, where n > 1) – Temperatures (K) at each time in tsteps
tsteps (NumPy array of floats with length n, where n > 1) – Array of times (yrs BP) in chronological (descending) order. First time is start of first timestep, last time is end of last timestep. Each pair of adjacent times composes a timestep
dpar (float) – Etch figure length (micrometers)
constants (dict) – Dictionary of constants associated with annealing model being used (default: KETCHAM_99_FC)
get_lengths (bool) – Boolean indicating whether or not to also calculate distribution of c-axis projected fission track lengths within an apatite that experienced this time-temperature history
make_graph (bool) – Boolean indicating whether or not to return length and frequency series for plotting length distribution (default: False). Only relevant if get_lengths is True.
l_num (int) – Nonnegative integer specifying how many lengths to use for the frequency series and internal calculations (default: 100). A value of at least 50 is recommended for calculating the standard deviation and mean of the mixed distribution; larger values produce smoother curves when plotting the frequency series. Only relevant if get_lengths is True

Returns:

age (float) – Model age (Myrs BP) of apatite that experienced the given time-temperature history. If get_lengths is False, only age is returned.
length_results (optional) (tuple) – Tuple describing the length distribution produced by calc_l_dist. Only included if get_lengths is True. Contains (in order) the mean (float) and standard deviation (float) of the distribution. If make_graph is also true, includes lengths (NumPy array of floats) and frequencies at which each length occurs in the distribution (NumPy array of floats)

gdtchron.aft.g(r, constants={'alpha': -0.12327, 'beta': -11.988, 'c0': -19.844, 'c1': 0.38951, 'c2': -51.253, 'c3': -7.6423, 'l_intercept': 15.72, 'l_slope': 0.35, 'r_kappa_sum': 1.0})[source]

Implement the length transform from Ketcham et al. (1999) (Equation 6).

Note: g is undefined for reduced lengths of 0 and may throw overflow exceptions for reduced lengths that are very close to 0. With the Ketcham et al. (1999) constants, this exception can occur for reduced lengths below 0.0007.

Parameters:

r (float or NumPy array of floats) – Reduced length (unitless)
constants (dict) – Dictionary of constants associated with annealing model being used (default: KETCHAM_99_FC)

Returns:

Result of length transform

Return type:

float or NumPy array of floats

gdtchron.aft.get_equiv_time(r_initial, temperature, constants={'alpha': -0.12327, 'beta': -11.988, 'c0': -19.844, 'c1': 0.38951, 'c2': -51.253, 'c3': -7.6423, 'l_intercept': 15.72, 'l_slope': 0.35, 'r_kappa_sum': 1.0})[source]

Calculate time it takes to reach a reduced length at given temperature.

This function solves Equation 5 from Ketcham (2005) for t (time).

Parameters:

r_initial (NumPy array of floats) – Reduced length (unitless)
temperature (float) – Temperature (K)
constants (dict) – Dictionary of constants associated with annealing model being used (default: KETCHAM_99_FC)

Returns:

Time(s) (in seconds) it would take to anneal to the given reduced length(s) if temperature remained constant

Return type:

NumPy array of floats

gdtchron.aft.get_l_stdev(l_c)[source]

Get standard deviation of length distribution based on its mean length.

This function follows equation from Figure 6b in Ketcham (2005). Assumes lengths are normally distributed around mean.

Parameters:: l_c (float (or NumPy array of floats)) – Mean c-axis projected length(s) (micrometers)
Returns:: Standard deviation(s) (micrometers) corresponding to each l value
Return type:: float (or NumPy array of floats)

gdtchron.aft.get_next_r(temperature, time_annealed, constants={'alpha': -0.12327, 'beta': -11.988, 'c0': -19.844, 'c1': 0.38951, 'c2': -51.253, 'c3': -7.6423, 'l_intercept': 15.72, 'l_slope': 0.35, 'r_kappa_sum': 1.0})[source]

Calculate reduced lengths after annealing over a given time period.

This function solves Equation 5 from Ketcham (2005) for r (reduced length).

Parameters:

temperature (float) – Temperature (K)
time_annealed (NumPy array of floats) – How long the crystal annealed at a given temperature (s)
constants (dict) – Dictionary of constants associated with annealing model being used (default: KETCHAM_99_FC)

Returns:

Mean reduced length(s) of fission tracks that annealed at the given temperature for the given period(s) of time

Return type:

NumPy array of floats

gdtchron.aft.l_conversion(r, dpar, constants={'alpha': -0.12327, 'beta': -11.988, 'c0': -19.844, 'c1': 0.38951, 'c2': -51.253, 'c3': -7.6423, 'l_intercept': 15.72, 'l_slope': 0.35, 'r_kappa_sum': 1.0})[source]

Convert reduced length r to mean c-axis projected length l_c.

Parameters:

r (float (or NumPy array of floats)) – reduced length (unitless)
dpar (float) – Etch figure length (micrometers) for the apatite
constants (dict) – Dictionary of constants associated with annealing model being used (default: KETCHAM_99_FC)

Returns:

Mean c-axis projected length(s) (micrometers) corresponding to each r value

Return type:

float or NumPy array of floats

gdtchron.aft.r_to_rho(r)[source]

Convert final reduced lengths to fission track densities.

This function uses Equation 13 from Ketcham (2005).

Parameters:: r (NumPy array of floats) – Mean reduced lengths (unitless) of fission tracks for the specific apatite for each point on the apatite’s time-temperature path. These values should already be converted to account for dpar variations.
Returns:: Normalized fission track densities corresponding to each interval on the time-temperature path. Any reduced lengths below 0.13 can’t be observed, so for intervals with a mean reduced length below 0.13, their corresponding fission track densities are set to 0.
Return type:: NumPy array of floats