#!/usr/bin/env python3
"""Independently remeasure saved synthetic ring positions and photon counts.
Standard-library Python3.10+. No Lab imports, private checkout or network required.
This is a frozen-record analysis, not a re-generation of the many-body evolution.
"""
import argparse
import hashlib
import json
import math
from pathlib import Path

HERE = Path(__file__).resolve().parent
ROOT = HERE.parents[2] if len(HERE.parents)>=3 else HERE
REPOSITORY_INPUT = ROOT/'docs/astronomy/ac10/synthetic-results.json'
DEFAULT_INPUT = REPOSITORY_INPUT if REPOSITORY_INPUT.exists() else HERE/'synthetic-results.json'


def measure(entry):
    record = entry['record']
    width = .28
    area = .5*6.8*16/len(record['particlePositions'])
    maximum_difference = 0
    transmission_estimates = []
    for point in record['points']:
        total = 0
        for particle in record['particlePositions']:
            dx = particle['x']-point['x']
            dy = particle['y']-16*math.floor((particle['y']+8)/16)
            if abs(dx)<=5*width and abs(dy)<=5*width:
                total += math.exp(-(dx*dx+dy*dy)/(2*width*width))/(2*math.pi*width*width)
        depth = area*total
        maximum_difference = max(maximum_difference,abs(depth-point['normalDepth']))
        predicted = math.exp(-depth/math.sin(math.radians(record['settings']['openingAngleDegrees'])))
        if abs(predicted-point['transmission'])>1e-10:raise ValueError('Saved ideal transmission differs from independent positional column')
        estimate = (point['photons']-4)/400
        if abs(estimate-point['estimatedTransmission'])>1e-12:raise ValueError('Saved photon estimate mismatch')
        transmission_estimates.append(estimate)
    if maximum_difference>1e-10:raise ValueError('Independent density integral mismatch')
    return {'preset':entry['preset'],'seed':record['seed'],'snapshot_time_inverse_omega':record['snapshotTime'],'tracers':len(record['particlePositions']),'maximum_depth_recalculation_error':maximum_difference,'mean_transmission':sum(transmission_estimates)/len(transmission_estimates),'minimum_measured_transmission':min(transmission_estimates),'nonpositive_net_samples':sum(value<=0 for value in transmission_estimates)}


def main():
    parser = argparse.ArgumentParser(description=__doc__)
    parser.add_argument('--input',type=Path,default=DEFAULT_INPUT)
    parser.add_argument('--out',type=Path,default=ROOT/'docs/astronomy/ac10/synthetic-replay.json' if REPOSITORY_INPUT.exists() else HERE/'synthetic-replay-rerun.json')
    args = parser.parse_args()
    raw = args.input.read_bytes()
    data = json.loads(raw)
    result = {'kind':'independent reanalysis of saved synthetic positions and detector counts; not measured Cassini particles','input_sha256':hashlib.sha256(raw).hexdigest(),'parameters':{'kernel_width_scene_units':.28,'scene_unit_meters':100,'unocculted_photons_per_sample':400,'known_background_photons':4},'results':[measure(entry)for entry in data['results']]}
    args.out.parent.mkdir(parents=True,exist_ok=True)
    args.out.write_text(json.dumps(result,indent=2,allow_nan=False)+'\n')
    print(json.dumps(result,indent=2))


if __name__=='__main__':main()