"""Compute finer whole states, retain all birth intervals and exact numeric refinement summaries."""
import signal
signal.alarm(120)
import argparse,hashlib,json
from pathlib import Path
import numpy as np
from evolution import simulate
from io_helpers import PACKAGE,fixed_inputs,family,read_json,write_json
from compare import finite_rows,log_probabilities
from observation_quadrature import nodes_from_saved_masses,quadrature_logdensity

def main():
 parser=argparse.ArgumentParser();parser.add_argument('--records',type=Path,default=Path('results/records'));parser.add_argument('--output',type=Path,default=Path('results/reference'));parser.add_argument('--real-data',type=Path);args=parser.parse_args()
 for case in family():
  path=args.records/(case['name']+'.npz');coarse=read_json(path);fine=simulate(case,fixed_inputs(),.0025,True);a=np.array(coarse['states']);indices=np.rint(np.array(coarse['timeGyr'])/.0025).astype(int);b=np.array(fine['states'])[indices]
  metrics={'bulk':float(np.max(np.abs(a[:,:4]-b[:,:4]))),'elementScaled':float(np.max(np.abs(a[:,4:]-b[:,4:])/np.maximum(np.abs(b[:,4:]),1e-8))),'abundanceDex':float(np.max(np.abs(np.array(coarse['abundances'])-np.array(fine['abundances'])[indices]))),'perFieldAbsolute':np.max(np.abs(a-b),axis=0).tolist()}
  datasets=[('synthetic',PACKAGE/'inputs/own-observations.json')];
  if args.real_data:datasets.append(('Hillconditional',args.real_data))
  prepared4=nodes_from_saved_masses(fine,fixed_inputs(),4);prepared8=nodes_from_saved_masses(fine,fixed_inputs(),8);observation=[]
  for label,data_path in datasets:
   rows=finite_rows(read_json(data_path));data_sha=hashlib.sha256(data_path.read_bytes()).hexdigest()
   for rho in [-.7,0,.7]:
    for axis,shift in [('none',0),('feh',-.1),('feh',.1),('mgfe',-.1),('mgfe',.1)]:
     fshift=shift if axis=='feh'else 0;mshift=shift if axis=='mgfe'else 0
     pc=log_probabilities(coarse,rows,rho,fshift,mshift);pf=log_probabilities(fine,rows,rho,fshift,mshift);p4,total4=quadrature_logdensity(fine,fixed_inputs(),rows,rho,4,fshift,mshift,prepared4);p8,total8=quadrature_logdensity(fine,fixed_inputs(),rows,rho,8,fshift,mshift,prepared8)
     observation.append({'data':label,'dataSha256':data_sha,'rho':rho,'shiftAxis':axis,'shiftDex':shift,'adjacentBirthBinMeanAbsoluteLogDensity':float(np.mean(np.abs(pc-pf))),'productionToPositiveGL8MeanAbsoluteLogDensity':float(np.mean(np.abs(pc-p8))),'GL4to8MeanAbsoluteLogDensity':float(np.mean(np.abs(p4-p8))),'analyticFormedNormalizationDifference':abs(total8-case['formedMass'])})
  # Retain the full fine cohort distribution and a separate actualfine-state hash.
  result={key:fine[key]for key in ['schema','name','parameters','saveEveryGyr','cohortFields','cohorts','minimumGas','units']};result.update(schema='chemical-reader-reference-cohorts-v1',fullFineRecordSha256=hashlib.sha256(json.dumps(fine,separators=(',',':'),allow_nan=False).encode()).hexdigest(),coarseRecordSha256=hashlib.sha256(path.read_bytes()).hexdigest(),refinement=metrics,observationRefinement=observation,observationReference='PositiveGL4/8integrals of analyticbirthrate over piecewise-linear finegas/elementmasses, before logarithms. Independentobservationoperator; same chemicalsolution.',earlyPositiveSupport={'minimumSavedGas':fine['minimumGas'],'firstFeFraction':fine['states'][0][5]/fine['states'][0][0],'firstMgFraction':fine['states'][0][4]/fine['states'][0][0]},fullFineBulkArraysArchived=False)
  write_json(args.output/(case['name']+'.npz'),result)
  if any(row['adjacentBirthBinMeanAbsoluteLogDensity']>=.01 or row['productionToPositiveGL8MeanAbsoluteLogDensity']>=.01 or row['GL4to8MeanAbsoluteLogDensity']>=.01 for row in observation):raise ValueError('Observation.01meanlogdensity gate failed; retainedactualrefinement')
  if metrics['bulk']>1e-6 or metrics['elementScaled']>.001 or metrics['abundanceDex']>.001:raise ValueError('Whole-state refinement gate failed')
 print('PASS:four finer whole states computed;cohorts and honestfieldmaxima retained')
if __name__=='__main__':main()