#!/usr/bin/env python3
"""Truth in Stone (2026), released under CC0 1.0.
Was the Bimini Road site dry land 10,000-12,000 years ago?

Data: Peltier & Fairbanks (2006) Barbados coral sea-level record, NOAA NCEI
https://www.ncei.noaa.gov/pub/data/paleo/coral/atlantic/barbados2006-noaa.txt
(dataset doi:10.25921/wtt6-rz64; paper doi:10.1016/j.quascirev.2006.04.010).

Acropora palmata grows in the top ~5 m of water, so for each A. palmata sample
sea level lay between its uplift-corrected depth and ~5 m above it. We take the
shallowest plausible sea level (depth + 5 m) as a conservative upper bound and
compare it with the depth of the Bimini Road blocks (4.6-5.5 m below today).
Assumes Bimini and Barbados saw similar relative sea level once Barbados uplift
is removed; glacial-isostatic differences across the Caribbean are a few metres.
"""
import urllib.request, sys
URL = "https://www.ncei.noaa.gov/pub/data/paleo/coral/atlantic/barbados2006-noaa.txt"
HABITAT = 5.0          # m, upper growth range of A. palmata
ROAD = (-5.5, -4.6)    # m, depth of the Bimini Road blocks

try:
    lines = open(sys.argv[1]).read().splitlines() if len(sys.argv) > 1 else \
        urllib.request.urlopen(URL, timeout=60).read().decode().splitlines()
except Exception as e:
    sys.exit(f"could not read data: {e}")
rows = [l.split("\t") for l in lines if l and not l.startswith("#")][1:]
pal = [(float(r[5]), float(r[4])) for r in rows if r[5].strip() and "palmat" in r[2]]

for lo, hi in [(8500, 9500), (9500, 10500), (10500, 11500), (11500, 12500)]:
    sel = [s for a, s in pal if lo <= a < hi]
    if not sel:
        continue
    top = max(sel) + HABITAT
    print(f"{lo/1000:.1f}-{hi/1000:.1f} ka: n={len(sel)}, corrected coral depths {min(sel):.1f} to {max(sel):.1f} m; "
          f"sea level no higher than about {top:.1f} m -> Road top ({ROAD[1]} m) was {ROAD[1]-top:.0f} m or more above sea level")
