"""Independent planning reference. Requires astropy; never downloads time tables.

Download targets.json and nasa-ps-source.json beside this file. Run explicitly:
python reproduce.py --target planet-hat-p-1-b --utc 2026-10-25T00:00:00
The UTC conversion uses Astropy's builtin Earth ephemeris and its available
leap-second table, independently of the browser's Astronomy Engine.
"""
import argparse, hashlib, json, math
from pathlib import Path

def predict(target, epoch, sigma_new=60):
    a=target['sigmaEpochDays']*86400;b=abs(epoch)*target['sigmaPeriodDays']*86400
    prior=math.hypot(a,b)
    return {'epoch':epoch,'predictedBjdTdb':target['epochBjdTdb']+epoch*target['periodDays'],
      'sigmaRhoZeroSeconds':prior,'sigmaCovarianceUpperBoundSeconds':a+b,
      'conditionalPosteriorSigmaSeconds':prior*sigma_new/math.hypot(prior,sigma_new),
      'conditionalInformationGainBits':.5*math.log2(1+prior*prior/(sigma_new*sigma_new))}

def astropy_bjd(utc,target,latitude=41.0082,longitude=28.9784,height=40):
    from astropy.time import Time
    from astropy.coordinates import SkyCoord, EarthLocation, solar_system_ephemeris
    from astropy.utils import iers
    import astropy.units as u
    iers.conf.auto_download=False
    site=EarthLocation.from_geodetic(longitude*u.deg,latitude*u.deg,height*u.m)
    time=Time(utc,scale='utc',location=site)
    coord=SkyCoord(ra=target['raDeg']*u.deg,dec=target['decDeg']*u.deg,frame='icrs')
    with solar_system_ephemeris.set('builtin'):
        return float((time.tdb+time.light_travel_time(coord,kind='barycentric')).jd)

def main():
    p=argparse.ArgumentParser();p.add_argument('--catalogue',type=Path,default=Path(__file__).with_name('targets.json'))
    p.add_argument('--target',default='planet-hat-p-1-b');p.add_argument('--utc',default='2026-10-25T00:00:00')
    p.add_argument('--latitude',type=float,default=41.0082);p.add_argument('--longitude',type=float,default=28.9784)
    p.add_argument('--height',type=float,default=40);p.add_argument('--measurement-sigma',type=float,default=60)
    a=p.parse_args();d=json.loads(a.catalogue.read_text());raw=a.catalogue.with_name('nasa-ps-source.json').read_bytes()
    assert hashlib.sha256(raw).hexdigest()==d['source']['sha256']
    row=next(x for x in d['targets']if x['id']==a.target)
    bjd=astropy_bjd(a.utc,row,a.latitude,a.longitude,a.height)
    n=math.ceil((bjd-row['epochBjdTdb'])/row['periodDays'])
    result=predict(row,n,a.measurement_sigma)
    print(json.dumps({'targetId':row['id'],'timeInputUTC':a.utc,'independentAstropyBjdTdb':bjd,
      'nextEpochPrediction':result,'sourceSha256':d['source']['sha256'],
      'limits':['Linear ephemeris, source duration fixed','Ranking assumes unpublished covariance equals zero','No real observation or instrument request','Future time depends on available leap-second table']},indent=2))

if __name__=='__main__':main()