#!/usr/bin/env python3
"""Recompute historical fixed-window approximations, not provider reliability."""
import argparse
import csv
import hashlib
import io
import json
import math
from decimal import Decimal, localcontext
from pathlib import Path

D = Decimal
YEAR_HOURS = D(8760)
NOTEBOOK_COMMIT = "37f758d64240cfb42a6d6fd41ff4ca16b8499bb3"


def annualize(probability, periods):
    if not D(0) <= probability <= D(1) or periods <= 0:
        raise ValueError("Probability must be in [0, 1]; periods must be positive")
    if probability in (D(0), D(1)):
        return probability
    return D(1) - ((D(1) - probability).ln() * periods).exp()


def compare(rate, hours):
    rate, hours = D(str(rate)), D(str(hours))
    if not rate.is_finite() or not hours.is_finite():
        raise ValueError("Inputs must be finite")
    if not D(0) <= rate <= D(1) or not D(0) < hours <= YEAR_HOURS:
        raise ValueError("Rate must be 0..1/year and window 0..8760 hours")
    with localcontext() as ctx:
        ctx.prec = 70
        periods = YEAR_HOURS / hours
        mu = D(20) * rate / periods
        exact_four = (-mu).exp() * mu**4 / D(24)
        # Positive tail summation avoids subtracting a near-one CDF.
        tail, term = exact_four, exact_four
        k = 4
        while term:
            k += 1
            term *= mu / D(k)
            previous = tail
            tail += term
            if tail == previous:
                break
        p = D(1) - (-rate / periods).exp()
        distinct = sum((D(math.comb(20, k)) * p**k * (1-p)**(20-k)
                        for k in range(4, 21)), D(0))
        return {
            "rate_per_drive_year": rate, "window_hours": hours,
            "periods_per_365_day_year": periods,
            "poisson_exact_four_per_window": exact_four,
            "poisson_at_least_four_per_window": tail,
            "binomial_distinct_at_least_four_per_window": distinct,
            "poisson_exact_four_annual": annualize(exact_four, periods),
            "poisson_at_least_four_annual": annualize(tail, periods),
            "binomial_distinct_at_least_four_annual": annualize(distinct, periods),
            "printed_periods": D(56) if hours == 156 else D(28) if hours == 312 else periods,
            "printed_poisson_annual": annualize(exact_four, D(56) if hours == 156 else D(28) if hours == 312 else periods),
        }


def build():
    scenarios = [("historical-inputs", "0.0041", 156),
                 ("without-recovery-credit", "0.0081", 156),
                 ("illustrative-two-percent", "0.02", 156),
                 ("double-window", "0.0041", 312)]
    rows = [{"id": name, **{key: str(value) for key, value in compare(rate, hours).items()}}
            for name, rate, hours in scenarios]
    with localcontext() as ctx:
        ctx.prec = 70
        p, count, years = D("1e-11"), D(10000000), D(10000)
        expected = count * p
        independent_any = annualize(p, count)
        ten_thousand_year_any = annualize(p, count * years)
    return {
        "schemaVersion": 1,
        "kind": "deterministic-model-arithmetic",
        "generatedOn": "2026-10-15",
        "generatorSha256": hashlib.sha256(Path(__file__).read_bytes()).hexdigest(),
        "precisionDecimalDigits": 70,
        "notebookCommit": NOTEBOOK_COMMIT,
        "assumptions": [
            "20 shards; recovery requires 17; constant independent per-shard hazard.",
            "Fixed windows are treated as independent, with no repair dynamics across window boundaries.",
            "Poisson counts events; binomial counts distinct shards failing at least once.",
            "Fractional 8760/hours exponent is a window-model convention, not a simulated partial window.",
            "No observed loss sample, confidence interval, S3 implementation or provider certification."
        ],
        "rows": rows,
        "expectationExample": {
            "objects": str(count), "annual_probability_per_object": str(p),
            "expected_losses_in_one_year": str(expected),
            "probability_any_loss_one_year_if_independent": str(independent_any),
            "probability_any_loss_10000_years_if_independent": str(ten_thousand_year_any),
            "note": "Same object-year independence assumption throughout; extrapolation illustrates expectation, not a forecast."
        }
    }


def main():
    parser = argparse.ArgumentParser()
    parser.add_argument("--output-dir", type=Path, required=True)
    args = parser.parse_args()
    data = build()
    args.output_dir.mkdir(parents=True, exist_ok=True)
    (args.output_dir / "durability-model.json").write_text(json.dumps(data, indent=2) + "\n")
    output = io.StringIO(newline="")
    writer = csv.DictWriter(output, fieldnames=list(data["rows"][0]), lineterminator="\n")
    writer.writeheader()
    writer.writerows(data["rows"])
    (args.output_dir / "durability-model.csv").write_text(output.getvalue())


if __name__ == "__main__":
    main()