Greetings community,
I am developing a age-depth model for radiocarbon dating. And there is a need to ‘calibrate’ radiocarbon determination (BP) into Calibrated date (cal AD). I have a .csv file containing information about calibration curve (first two rows given below).
CAL BP |
14C age |
Sigma |
55000 |
50100 |
1024 |
54980 |
50081 |
1018 |
The radiocarbon determination has 14C age +/ sigma (can be seen on y-axis in image). I want to achieve graph of the x-axis, given the calibration curve csv and radiocarbon determination.
Thank you.
What have you tried already?
Greetings Ricardo,
I have tried to implement it using scipy.stats.
# converts 14C age to calendar age
def calib(c14_age, c14_error, calibration_curve, cutoff=0.01):
probs = np.column_stack((
calibration_curve[:, 0], # calendar ages
norm.pdf(calibration_curve[:, 1], c14_age,
np.sqrt(calibration_curve[:, 2]**2 + c14_error**2))
))
probs[probs[:, 1] < 0, 1] = 0
# normalization
total = np.sum(probs[:, 1])
if total > 0:
probs[:, 1] = probs[:, 1] / total
# apply cutoff
above = np.where(probs[:, 1] / np.max(probs[:, 1]) > cutoff)[0]
if len(above) > 5:
filtered_probs = probs[np.min(above):np.max(above) + 1].copy()
filtered_probs[filtered_probs[:, 1] < 0, 1] = 0
if np.sum(filtered_probs[:, 1]) > 0:
filtered_probs[:, 1] = filtered_probs[:, 1] / np.sum(filtered_probs[:, 1])
probs = filtered_probs
return probs