"""Reproduce reference probabilities and transparent sensitivity scenarios. Standard library only."""
from pathlib import Path
import json
import math

HERE = Path(__file__).resolve().parent
data = json.loads((HERE / 'inputs.json').read_text())

def at_least_one(p, n):
    return -math.expm1(n * math.log1p(-p))

rows = []
for n, published in zip(data['years'], data['published_percentages']):
    for p, reference in zip(data['annual_probabilities'], published):
        calculated = at_least_one(p, n) * 100
        assert abs(calculated - reference) <= 0.5
        rows.append({'years': n, 'annual_probability': p, 'reference_percent': reference,
                     'calculated_percent': calculated, 'difference_percentage_points': calculated-reference})
p, n = 0.01, 30
two_or_more = 1 - (1-p)**n - n*p*(1-p)**(n-1)
scenario = data['sensitivity_scenario']
annual = [scenario['start_probability'] + i*(scenario['end_probability']-scenario['start_probability'])/(n-1) for i in range(n)]
changing = 1-math.prod(1-value for value in annual)
report = {'reference_comparison': rows,
          'one_percent_30_year_at_least_one_percent': at_least_one(p,n)*100,
          'one_percent_30_year_at_least_two_exceedance_years_percent': two_or_more*100,
          'one_percent_30_year_expected_exceedance_years': n*p,
          'illustrative_rising_probability_30_year_percent': changing*100,
          'constant_two_percent_30_year_percent': at_least_one(.02,n)*100,
          'interpretation': 'Reference consistency check plus analytic scenarios; not observed flood data or a local forecast.'}
(HERE/'results.json').write_text(json.dumps(report,indent=2)+'\n')
print(json.dumps({k:v for k,v in report.items() if k != 'reference_comparison'},indent=2))